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This  study  investigates  the  problems  of  model  refinement  and  structural  damage 
assessment.  The  essence  of  the  model  refinement  problem  is  to  adjust  finite  element  models 
(FEMs)  of  structures  with  the  intent  of  producing  a  correlation  between  experimental  and 
analytical  modal  properties.  Recently,  the  framework  of  model  refinement  has  been  adopted 
to  determine  the  location  and  extent  of  structural  damage.  Damage  will  result  in  changes  to 
the  modal  properties  of  the  healthy  structure.  A  further  refinement  of  an  accurate  FEM  of  the 
structure  using  damaged  modal  parameters  is  expected  to  generate  adjustments  to  the  FEM  at 
locations  associated  with  damage.  Analysis  of  these  adjustments  can  then  be  used  to  assess 
damage.  In  this  investigation,  four  algorithms  relevant  to  the  subjects  of  model  refinement 
and  damage  detection  are  presented.  The  development  of  a  model  refinement  algorithm  with 
its  basis  in  inverse  theory  is  first  presented.  The  algorithm,  termed  the  inverse/hybrid 
approach,  is  illustrated  in  a  comparative  computer  simulated  study.  Next,  an  existing 
eigenstructure  assignment  model  refinement  algorithm  is  improved  to  better  approach  the 
damage  assessment  problem.  The  enhanced  algorithm  is  evaluated  and  compared  to  other 
techniques  using  simulated  and  experimental  data.  The  algorithm  is  shown  to  perform  well 
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in  assessing  damage  and  refining  FEMs.  A  damage  location  algorithm  that  bypasses  the 
general  framework  of  model  refinement  is  discussed.  The  damage  location  algorithm, 
termed  the  subspace  rotation,  is  similar  to  the  modal  force  error  criteria  proposed  by  several 
researchers.  Greater  insight  to  the  modal  force  error  criteria,  along  with  a  new  viewpoint  that 
reduces  the  effects  of  measurement  noise,  is  discussed.  Furthermore,  an  efficient  damage 
extent  algorithm  based  on  a  minimum  rank  perturbation  theory  (MRPT)  is  developed.  The 
formulation  of  the  MRPT  is  consistent  with  the  effect  of  most  structural  damage  on  FEMs. 
The  characteristics  of  the  subspace  rotation  algorithm  and  the  minimum  rank  perturbation 
theory  are  illustrated  using  simulated  and  experimental  testbeds.  The  decomposition  of  the 
damage  assessment  problem  into  location  and  extent  subproblems  is  shown  to  be 
advantageous  from  both  for  computational  efficiency  and  for  engineering  insight. 
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CHAPTER  1 
INTRODUCTION 


1.1    Finite  Element  Model  Refinement 

1.1.1     Overview 

An  important  aspect  in  the  design  process  of  a  structure  is  the  evaluation  of  its 
performance  under  expected  dynamic  loading  conditions.  Dynamic  performance  can  be 
analyzed  by  using  either  analytical  or  experimental  techniques.  Experimental  analyses  are 
generally  very  costly,  time  consuming  and  can  encounter  technical  difficulties.  One  example 
of  a  cormnonly  encountered  technical  barrier  is  the  ground  base  reproduction  of  a  weightless 
environment  for  the  testing  of  space  structures.  This  shortcoming  and  the  advent  of 
computers  have  sparked  a  growing  interest  in  the  uses  of  analytical  techniques.  This  type  of 
analysis  utilizes  a  mathematical  model  of  the  actual  structure  for  computer  simulated 
evaluation  of  structural  performance.  Presently,  most  analytical  models  used  in  industry  are 
finite  element  models  generated  by  using  Finite  Element  Methods.  The  accuracy  of  a  finite 
element  model  (FEM)  can  be  improved  to  some  extent  by  increasing  the  number  of  degrees 
of  freedom  (DOFs)  included  in  the  model.  However,  accuracy  of  the  FEM  may  still  be 
lacking  due  to  uncertainties  in  material  properties,  geometric  layout  and  fabrication  induced 
errors.  Thus,  it  is  essential  to  "validate"  these  FEMs  prior  to  their  acceptance  as  a  basis  for 
analysis.  One  way  to  validate  them  is  to  compare  their  modal  properties  (eigenvalues/ 
eigenvectors)  with  the  measured  modal  parameters  of  the  actual  structure  as  obtained  from 
experimental  modal  analysis  (EM A).  A  FEM  is  acceptable  when  these  two  sets  of  modal 
properties  are  in  agreement.  Unfortunately,  this  agreement  rarely  occurs.  As  a  result,  the 
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FEM  must  be  adjusted  to  produce  a  correlation  of  analytical  and  experimental  modal 
properties.  An  overview  of  this  procedure  is  schematically  illustrated  in  Figure  1.1 
(Zimmerman  and  Smith,  1992). 
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Figure  1 . 1        Overview  of  Finite  Element  Model  Refinement 
(Zimmerman  and  Smith,  1992) 


In  the  past,  the  process  of  adjusting  FEMs  was  performed  on  an  ad-hoc  basis  aided  by 
engineering  experience.  This  practice  was  naturally  time  consuming  and  in  most  cases 
inadequate  for  large-order,  complex  FEMs  that  are  typically  needed  for  accurate  dynamic 
modelling.  In  the  past  years,  growing  interest  has  been  focused  on  developing  systematic 
procedures  to  produce  correlated  FEMs.  These  efforts  have  resulted  in  the  development  of  a 
large  number  of  algorithms.    Finite  element  model  adjustment  procedures  have  been 


commonly  referred  to  in  the  literature  as,  FEM  refinement,  FEM  adjustment,  FEM 
correction,  FEM  correlation,  and  FEM  identification. 

1.1.2     Literature  Survey 

Comprehensive  literature  surveys  covering  a  large  portion  of  the  work  that  addresses  the 
model  refinement  problem  can  be  found  in  the  book  chapter  by  Zimmerman  and  Smith 
(1992),  and  in  the  papers  by  Ibrahim  and  Saafan  (1987)  and  Heylen  and  Saas  (1987).  For 
completeness  and  to  properly  underline  the  objective  of  the  current  study,  it  is  appropriate  to 
include  a  brief  survey  of  the  development  made  in  this  area. 

The  concept  of  using  experimental  modal  data  in  analytical  studies  was  initiated  by  an 
early  work  presented  by  Rodden  (1967).  In  his  work,  Rodden  explored  the  possibilities  of 
generating  mass  and  stiffness  property  matrices  by  using  experimentally  measured  modal 
data.  The  experimentally  generated  mass  and  stiffness  matrices  were  nonsymmetric.  Brock 
(1968)  improved  the  work  of  Rodden  by  proposing  a  strategy  to  insure  symmetry  of  property 
matrices. 

The  essence  of  the  model  refinement  concept,  as  adopted  by  most  researchers,  is  to 
modify  the  finite  element  model  (FEM)  of  systems  with  the  intention  of  producing  a 
correlation  of  test  and  analytical  modal  parameters.  In  general,  this  task  has  been  approached 
through  two  different  philosophies.  The  first  amounts  to  modifying  globally  the  finite 
element  property  matrices  of  the  system.  This  approach  is  commonly  referred  to  as  the 
property  matrix  update.  The  alternative  approach  is  to  individually  correct  parameters  of 
each  element  of  the  finite  element  model.  This  will  be  addressed  here  as  the  physical 
parameter  update. 

Property  matrix  update  algorithms  can  be  classified  into  two  subclasses  of 
methodologies:  optimal-matrix  update  and  control-based  eigenstructure  assignment 
techniques. 
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The  basic  philosophy  of  the  optimal-matrix  update  is  to  minimize  the  correction  to  the 
FEM  property  matrices  to  accomplish  the  analytical/test  modal  correlation.  The  pioneering 
work  in  this  area  can  be  credited  to  Baruch  and  Bar  Itzhack  (1978).  In  their  formulation,  by 
assuming  that  the  mass  matrix  is  correct,  the  refinement  of  the  stiffness  matrix  of  an 
undamped  FEM  is  cast  as  a  constrained  minimization  problem.  The  objective  of  their 
formulation  is  to  determine  the  minimal  Frobenius-norm  symmetric  stiffness  adjustment 
that  satisfies  the  eigenvalue  problem  in  terms  of  the  experimentally  measured  eigenvalues 
and  eigenvectors.  A  computationally  efficient  closed  form  solution  was  developed  for  the 
updated  stiffness  matrix.  Berman  and  Nagy  (1983)  extended  the  Baruch  and  Bar  Itzhack 
approach  to  the  refinement  of  both  the  mass  and  stiffness  matrices.  This  same  methodology 
was  further  extended  by  Fuh  et  al.  ( 1 984)  to  update  the  mass,  damping  and  stiffness  matrices 
of  damped  FEMs.  Fuh  and  his  colleagues  used  cross-orthogonality  relationships  to  correct 
the  mass  and  damping  matrices  and  a  constrained  minimization  problem  similar  to  the  one 
proposed  by  Baruch  and  Bar  Itzhack  to  determine  the  refined  stiffness  matrix.  The  problem 
of  adjusting  the  mass,  damping  and  stiffness  matrices  was  also  attempted  by  Hanagud  et  al. 
(1984).  In  their  approach,  all  three  property  matrices  of  a  nonproportionally  damped  FEM 
are  incorporated  in  the  constrained  minimization  problem. 

The  previously  stated  algorithms  do  not  preserve  the  load  path  (sparsity  pattern)  of  the 
original  analytical  property  matrices.  Kabe  (1985)  proposed  a  reformulafion  of  the  Baruch 
and  Bar  Itzhack  algorithm  that  constrains  the  updated  stiffness  matrix  to  preserve  the  load 
path  of  the  original  stiffness  matrix.  In  addition,  in  his  reformulation,  he  adopted  a  percent 
weighting  on  the  stiffness  changes  instead  of  mass  matrix  weighting  as  used  in  the  Baruch 
and  Bar  Itzhack  algorithm.  Kammer  ( 1 987)  proposed  an  alternative  solution  approach  to  the 
problem  defined  by  Kabe  that  uses  projection  matrix  theory  and  the  Moore-Penrose 
generalized  inverse.  The  algorithms  formulated  by  Kabe  and  Kammer  involve  an  excessive 
amount  of  computational  effort.  Another  alternative  and  more  efficient  formulation  of  the 


Kabe  problem  that  utilized  a  generalization  of  Marwill-Toint  updates  was  developed  by 
Smith  and  Beattie  (1991). 

The  other  subclass  of  property  matrix  updates  is  based  on  the  general  framework  of 
control-based  eigenstructure  assignment  algorithms.  The  essence  of  this  approach  is  to 
determine  pseudo-controllers  that  would  assign  the  experimentally  measured  modal 
parameters  to  the  original  analytical  FEM.  The  pseudo-controllers  are  then  translated  into 
matrix  adjustment  to  the  original  property  matrices.  This  approach  was  initially  formulated 
by  Inman  and  Minas  (1990)  to  adjust  the  damping  and  stiffness  matrices  of  the  original  FEM. 
In  their  formulation,  the  symmetry  of  the  adjusted  property  matrices  was  enforced  through 
an  iterative  process  that  involves  a  numerical  nonlinear  optimization  process.  Zimmerman 
and  Widengren  (1989,  1990)  proposed  a  formulation  that  replaces  the  iterative  process  of 
Inman  and  Minas  with  a  generalized  algebraic  Riccati  equation.  Further  details  about  this 
approach  are  treated  in  Chapter  4. 

The  alternative  philosophy  to  the  model  refinement  problem  is  the  physical  parameter 
update.  The  basic  idea  of  this  alternative  approach  consists  of  adjusting  some  or  all  structural 
physical  parameters  to  produce  test/analysis  modal  correlation.  Structural  parameters  are 
design  variables  such  as  material  densities,  modulus  of  elasticity,  cross-sectional  areas, 
element  lengths  and  others.  This  type  of  procedure  results  in  corrected  FEMs  that  are 
consistent  with  the  framework  of  the  original  FEM. 

Most  methodologies  that  have  adopted  this  philosophy  have  used  the  sensitivity 
derivatives  of  the  system  eigenvalues  and/or  eigenvectors  with  respect  to  the  structural 
parameters.  Generally,  the  refinement  process  amounts  to  solving  for  the  corrected 
structural  parameters  through  an  iterative  optimization  problem  directed  by  the  sensitivity 
derivatives.  Some  notable  work  in  this  area  was  accomplished  by  Collins  et  al.  (1974),  Chen 
and  Garba  (1980),  Adelman  and  Haftka  (1986),  Creamer  and  Hendricks  (1987),  Flanigan 
(1991),  Martinez  et  al.  (1991),  to  name  only  a  few.  In  these  studies,  different  sensitivity 
formulations  and  iterative  optimization  strategies  are  proposed.    A  direct  approach  that 


bypasses  the  use  of  closed-form  sensitivity  derivatives  was  presented  in  the  papers  of  Hajela 
and  Soeiro  (1990)  and  Soeiro  (1990). 

White  and  Maytum  (1976)  set  forth  another  alternative  methodology  that  can  be 
classified  in  the  physical  parameter  update  group.  In  their  approach,  adjustments  to  the 
stiffness  matrix  are  viewed  as  a  collection  of  known  submatrices,  elements  or  group  of 
elements.  Each  submatrix  is  multiplied  by  an  unknown  scaling  factor.  Correction  of  the 
original  stiffness  matrix  amounts  to  determining  the  scaling  factors  that  would  improve  the 
correlation  of  analytical  and  experimental  modal  properties.  The  scaling  factors  can  be 
physically  viewed  as  functions  of  the  stiffness  structural  parameters.  A  large 
underdetermined  system  of  equations  in  terms  of  the  scaling  factor  is  generated.  The  scaling 
factors  are  then  computed  by  using  a  pseudo-inverse.  An  improvement  to  the  White  and 
Maytum  approach  is  presented  in  the  paper  by  Lim  (1990). 

1.2    Structural  Damage  Assessment 

Structures,  in  general,  are  prone  to  structural  damage  during  their  service  lives,  that 
could  affect  their  overall  performance  and  could  result  in  catastrophic  failures.  This  is  of 
critical  concern  when  expensive  systems  and/or  human  lives  are  at  stake.  On  this  basis,  it  is 
highly  desirable  to  monitor  the  structural  health  of  systems  such  as  airplanes,  space  systems, 
bridges,  buildings,  oil  platforms,  among  others,  to  prevent  such  catastrophic  events. 
Currently,  the  most  common  structural  health  monitoring  methods  involve  visual 
inspections  supplemented  with  X-rays,  acoustic  emission,  magnetic  resonance  and 
ultrasonic  testing.  These  approaches  can  be  time  consuming,  costly  and  difficult  to  perform 
on  inaccessible  structural  components.  Additionally,  none  of  these  approaches  provide  a 
quantitative  assessment  of  the  magnitude  of  the  damage. 

In  response  to  these  shortcomings,  structural  health  monitoring  approaches  based  on 
the  framework  of  existing  model  refinement  techniques  have  been  recently  proposed  by 
several  researchers  (Smith  and  Hendricks,  1987;  Chen  and  Garba,  1988;  Ricles  and 


Kosmatka,  1992;  Zimmerman  and  Kaouk,  1992a,b;  Kaouk  and  Zimmerman,  1993a,b). 
These  techniques  utilize  the  vibration  signature  of  the  pre-  and  post-  damaged  structure,  in 
conjunction  with  an  analytical  model  of  the  original  structure,  to  assess  both  the  location  and 
extent  of  the  structural  damage.  The  pre-damage  modal  parameters  are  used  to  correct 
(refine)  the  original  finite  element  model  (FEM)  to  determine  an  "accurate"  reference 
baseline.  Once  damaged,  the  post-damage  modal  properties  of  the  structure  are  used  to 
further  "refine"  the  refined  analytical  model.  This  results  in  perturbations  to  the  refined 
analytical  model.  Analysis  of  the  perturbations  could  indicate  the  damage  location  and 
extent.  An  overview  of  the  application  of  model  refinement  algorithms  in  assessing 
structural  damage  is  summarized  in  Figure  1.2. 
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Figure  1 .2        Overview  of  Finite  Element  Model  Refinement 
Process  Used  for  Damage  Assessment 
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Notable  exceptions  to  the  direct  use  of  FEM  refinement  algorithms  to  the  damage 
detection  problem  are  the  work  of  Lin  (1990),  Ojalvo  and  Pilon  (1988),  and  Gysin  (1990).  In 
the  work  of  Lin,  a  flexibility  matrix  is  determined  using  experimental  data.  This  matrix  is 
then  multiplied  by  the  original  stiffness  matrix,  with  those  rows  and/or  columns  that  differ 
significantly  from  a  row  and/or  column  of  the  identity  matrix  indicating  which  degrees  of 
freedom  have  been  most  affected  by  the  damage.  It  is  then  assumed  that  damage  has  occurred 
in  structural  elements  connecting  those  degrees  of  freedom.  An  overview  of  the  work  by 
Ojalvo  and  Pilon  (1988),  and  Gysin  (1990)  is  given  in  Chapter  5. 

Although  the  problems  of  damage  assessment  and  model  refinement  are  related  to  one 
another,  they  have  quite  different  characteristics.  In  essence,  the  model  refinement  concept 
is  based  on  the  philosophy  of  the  minimum  change  to  the  original  FEM.  Naturally,  this 
minimum  change  constraint  has  a  tendency  of  smearing  the  changes  throughout  the  entire 
FEM.  This  aspect  is  inconsistent  with  the  effects  of  structural  damage  on  FEMs  which  are 
usually  localized  perturbations  of  possibly  large  magnitude.  In  fact,  structural  damage  often 
occurs  at  discrete  locations  and  only  affects  a  few  elements  of  the  FEM. 

1.3    Objective  of  the  Present  Study 

The  present  study  investigates  the  development  of  new  and  promising  model  refinement 
and  damage  detection  methodologies.  Although  considerable  research  has  been  done  in 
these  areas,  no  methodology  has  been  fully  successful  in  dealing  with  the  refinement  or  the 
damage  detection  problem  of  "real  life"  systems.  The  main  objective  of  this  investigation  is 
to  formulate  efficient  model  refinement  algorithms  that  are  consistent  with  the  effect  of 
structural  damage. 

In  chapter  2,  practical  concepts  and  issues  related  to  the  general  areas  of  model 
refinement  are  presented.  First,  the  concept  of  the  eigenvalue  problem  of  discrete  structural 
models  is  reviewed  with  emphasis  on  the  associated  cross-orthogonality  conditions.  A  brief 
discussion  of  modal  analysis  follows.  The  problem  associated  with  incomplete  eigenvector 
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measurements  is  then  investigated  and  already  existing  techniques  to  deal  with  this  problem 
are  presented.  Furthermore,  eigenvector  orthogonalization  techniques  useful  to  a  large 
number  of  model  refinement  algorithms  are  reported.  Finally,  an  iterative  approach  to  the 
problem  of  load  path  preservation  encountered  in  a  large  number  of  property  matrix  update 
algorithms  is  discussed. 

In  chapter  3,  the  author  proposes  and  formulates  a  new  model  refinement  algorithm. 
The  newly  developed  algorithm,  named  the  inverse/hybrid  method,  is  naturally  based  on  the 
inverse  problem.  The  model  type  under  consideration  is  proportionally  damped  and  the 
inaccuracies  of  the  analytical  model  is  assumed  to  be  due  to  errors  in  the  modelling  of  the 
damping  and  stiffness  properties.  A  complete  hybrid  set  of  modal  data  is  generated  by 
approximating  the  unmeasured  modal  information  by  the  corresponding  analytical  modes. 
Symmetry  of  the  corrected  stiffness  and  damping  matrices  is  enforced  by  mass 
orthogonalization  the  complete  hybrid  eigenvector  matrix.  An  orthogonalization  procedure 
that  assigns  more  confidence  on  the  measured  eigenvectors  is  proposed.  A  comparative 
study  of  the  inverse/hybrid  method  and  the  algorithm  proposed  by  Baruch  and  Bar  Itzhack 
(1978)  shows  that  both  algorithms  give  similar  results.  However,  it  is  shown  that  the 
inverse/hybrid  approach  is  less  computationally  efficient. 

Chapter  4  proposes  improvements  to  the  symmetric  eigenstructure  assignment  model 
refinement  algorithm  (SEAMRA)  formulated  by  Zimmerman  and  Widengren  (1989).  The 
author  develops  a  technique  to  enhance  eigenvector  assignability.  The  technique,  termed  the 
subspace  rotation  method,  is  based  on  rotating  the  achievable  eigenvector  subspaces  into  the 
experimental  eigenvectors.  The  subspace  rotation  method  results  in  both  a  decrease  in  the 
computational  burden  and  an  increase  in  the  accuracy  of  the  assigned  eigenvector.  The 
enhanced  SEAMRA  is  then  evaluated  and  compared  to  other  algorithms  using  both 
computer  simulated  and  experimental  testbeds.  It  is  shown  that  the  enhanced  algorithm  is 
suitable  for  damage  detection  applications. 
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In  Chapter  5,  an  efficient  damage  location  algorithm  that  utilizes  modal  data 
information  but  bypasses  the  general  framework  of  the  model  refinement  problem  is 
presented.  This  location  algorithm  is  an  outgrowth  of  the  subspace  rotation  method  used  to 
enhance  eigenvector  assignability  of  the  SEAMRA.  The  proposed  location  algorithm  is 
similar  to  the  modal  force  error  criteria  presented  by  several  researchers  (Ojalvo  and  Pilon, 
1988;Gysin,  1990).  Further  interpretation  ofthe  algorithm  operation  is  given.  Additionally, 
the  author  proposes  and  formulates  a  new  viewpoint  that  reduces  the  effect  of  measurement 
noise  for  certain  types  of  structures.  Based  on  this  formulation,  an  eigenvector  filtering 
algorithm  is  also  developed. 

Chapter  6  presents  the  formulation  of  computationally  attractive  damage  extent 
algorithms  that  are  based  on  a  minimum  rank  perturbation  theory  (MRPT)  developed  by  the 
author.  The  formulation  of  the  MRPT  is  consistent  with  the  effect  of  many  classes  of 
structural  damage  on  FEMs.  Several  MRPT  based  algorithms  are  formulated  to 
accommodate  structures  with  undamped,  proportionally  damped,  and  nonproportionally 
damped  FEMs.  For  each  type  of  FEM,  several  damage  scenarios  are  considered. 
Discussions  of  the  characteristics  and  properties  of  the  developed  algorithms  are  presented 
along  with  practical  issues  that  can  be  used  to  improve  their  performance. 

In  Chapter  7,  the  algorithms  developed  in  Chapters  5  and  6  are  demonstrated  and 
evaluated  using  both  computer  simulated  and  actual  experimental  data.  The  main  objective 
of  these  example  problems  is  to  illustrate  the  potential  of  these  algorithms  in  assessing 
structural  damage.  In  all  examples,  the  location  of  damage  is  first  determined  by  using  the 
algorithm  presented  in  Chapter  5.  An  MRPT  based  algorithm  is  then  utilized  to  assess  the 
extent  ofthe  damage.  Several  key  points  made  throughout  the  formulation  in  Chapters  5  and 
6  are  emphasized.  In  particular,  it  is  shown  that  the  damage  extent  calculations  can  be  greatly 
enhanced  by  making  use  of  the  damage  location  algorithm. 

In  Chapter  8,  summaries  and  conclusions  of  the  issues  discussed  in  this  study  are 
presented  along  with  some  suggestions  for  future  work. 


CHAPTER  2 

MATHEMATICAL  PRELIMINARffiS  AND  PRACTICAL  ISSUES  RELATED  TO 

THE  PROBLEMS  OF  MODEL  REFINEMENT  AND  DAMAGE  DETECTION 


2.1    Introduction 

In  this  chapter,  general  mathematical  preliminaries  and  practical  issues  relevant  to  the 
areas  of  model  refinement  and  damage  detection  are  discussed.  In  Section  2.2,  an  overview 
of  the  general  eigenvalue  problem  for  discrete  models  is  given.  Further,  the  concept  of  modal 
analysis  is  introduced  in  Section  2.3.  Section  2.4  deals  with  the  concept  of  incomplete 
degrees  of  freedom  measurement.  Two  alternative  approaches  are  discussed  as  possible 
solutions  to  the  incomplete  measurement  problem.  Two  eigenvector  orthogonalization 
techniques  are  discussed  in  Section  2.5.  Finally,  in  Section  2.6,  an  algorithm  to  preserve 
sparsity  in  updated  property  matrices  is  presented.  The  concepts  discussed  in  this  chapter 
will  be  frequently  referred  to  in  the  course  of  the  forthcoming  chapters. 

2.2   The  Eigenvalue  Problem  of  Discrete  Systems 

In  practice,  most  engineering  structures  are  continuous  systems  with  spatially 
distributed  material  properties.  The  vibration  motion  of  these  structures  in  terms  of  their 
distributed  properties  is  usually  governed  by  one  or  several  partial  differential  equations. 
With  complex  systems,  both  the  development  and  the  analysis  of  partial  differential 
equations  of  motion  are  tedious  and  in  many  cases  impossible.  These  shortcomings  and  the 
advent  of  digital  computers  have  motivated  the  development  of  approximate  modelling  of 
continuous  systems  in  order  to  simplify  the  equations  of  motion.  The  general  idea  behind 
these  approximations  is  to  represent  the  exact  distributed  model  of  a  system  by  a  discrete  one. 
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This  concept  is  known  as  spatial  discretization  which  eliminates  the  continuous  spatial 
dependence  of  the  distributed  properties.  The  discrete  model  is  of  finite  order  and  is 
described  by  a  finite  number  of  variables  known  as  degrees  of  freedom  (DOFs).  The  number 
of  DOFs  used  in  the  discrete  model  depends  on  the  desired  accuracy  in  representing  the 
continuous  model.  Commonly,  the  vibration  motion  of  systems  in  terms  of  their  discrete 
models  is  described  by  a  set  of  simultaneous  ordinary  differential  equations  that  are  usually 
simpler  to  develop  and  analyze  than  the  partial  differential  equations  of  the  continuous 
system.  One  of  the  most  commonly  used  approximate  discrete  modelling  techniques  is  the 
finite  element  method.  The  model  generated  by  a  finite  element  method  is  called  a  finite 
element  model  (FEM).  More  elaborate  discussions  of  the  concepts  of  continuous  and 
discrete  models  as  well  as  finite  element  methods  are  covered  in  detail  in  the  books  of 
Meirovitch  (1986,  1980),  Inman  (1989),  and  Hughes  (1987). 

Commonly,  the  free  vibration  motion  of  a  structure  in  terms  of  an  n^^  order  discrete 
model  is  represented  by  the  following  set  of  simultaneous  ordinary  differential  equations 

Mw(t)  +  Dw(t)  +  Kw(t)  =  0  (2.2.1) 

where  M,  D  and  K  are  termed,  respectively,  the  mass,  damping  and  stiffness  matrices.  They 
are  models  of  the  mass,  damping  and  stiffness  properties  of  the  structure.  Since  the  system 
model  is  order  n,  these  matrices  are  of  dimension  nxn  and  are  generally  real.  The  variable 
w(t)  represents  the  n  displacements  of  the  n-DOF  model  of  the  structure.  The  overdots 
represent  differentiation  with  respect  to  time.  The  mass  matrix,  M,  is  always  symmetric 
positive  definite.  The  stiffness  (K)  and  damping  (D)  matrices  of  nongyroscopic  and 
noncirculatory  systems  are  symmetric.  In  general,  the  modelling  of  the  mass  and  stiffness 
properties  of  the  structure  is  simpler  and  more  accurate  than  the  modelling  of  the  damping 
properties.  In  the  forthcoming  discussion,  it  is  assumed  that  the  system  under  consideration 
is  nongyroscopic  and  noncirculatory.  The  standard  solution  to  Eq.  (2.2.1)  is  of  the  form 


w(t)  =  ve^'  (2.2.2) 
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where  v  is  a  constant  nxl  vector  and  X,  is  a  constant  scalar.  Substituting  Eq.  (2.2.2)  into  Eq. 
(2.2.1)  and  dividing  the  resultant  equation  by  e    yields  the  condition 


Mv  X^  +  Dv  X  +  Kv  =  0 


(2.2.3) 


There  are  n  sets  of  nontrivial  complex  conjugate  solutions  (ki,  vj)  to  Eq.  (2.2.3).  Note  that 
since  the  property  matrices  (M,  D,  K)  are  real,  if  (A^,  Vj)  is  a  solution  set  to  Eq.  (2.2.3),  the 
complex  conjugate  of  that  set  is  also  a  solution.  The  problem  of  solving  for  these  solutions 
is  commonly  known  as  the  eigenvalue  problem  and  is  sometimes  referred  to  as  the 
characteristic  value  problem.  The  scalar  A,,  and  the  vector  v;  are  known,  respectively,  as  the 
eigenvalue  and  eigenvector  of  the  i^*^  mode  of  vibration  of  the  structure.  For  a  general 
damped  system,  eigenvalues  and  eigenvectors  are  both  complex.  Note  that  Eq.  (2.2.3)  can 
be  easily  rearranged  in  the  more  general  mathematical  eigenvalue  problem  format  as 


(2.2.4) 


where  Inxn  is  the  nxn  identity  matrix.  Equation  (2.2.4)  is  called  a  state  space  representation 
of  Eq.  (2.2.3). 

The  eigenvalue  and  eigenvector  can  be  related  to  some  physical  characteristic  properties 
of  structures.  For  that  matter,  the  i'*^  eigenvalue  is  written  as 


0 

^nxn 

■Xi" 

"Xi 

M-^K 

-  m-'d 

Xi^i 

=  ^i 

VjX 

Xi=    -^iODi+jtOi/T^ 


(2.2.5) 


where  j  =  /—  1  .  Note  that  in  writing  this  equation  it  is  assumed  that  the  system  is 
underdamped.  The  real  scalar  OO;  and  ^j  are  the  natural  frequency  and  damping  ratio  (or 
damping  factor),  respectively,  of  the  i"'  mode  of  the  structure.  The  eigenvector  Vj  indicates 
the  "shape"  of  vibration  of  the  i*  mode  of  the  structure.  The  sets  of  frequencies,  damping 
ratios  and  mode  shapes  are  sometimes  referred  to  as  modal  parameters. 
The  symmetric  nature  of  the  property  matrices  (M,  D,  K)  constrains  the  eigenvectors  to 
satisfy  some  cross-orthogonality  relationships.     For  the  purpose  of  discussing  these 
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cross-orthogonality    relationships,    consider    the    following    alternative    state    space 
representation  of  the  eigenvalue  problem  in  Eq.  (2.2.3) 


(2.2.6) 
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-  K 

VA  VA 

V      V 

A   0 

0   A 
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K 
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V      V 
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where         V  =      v,,  Vj,   .  . 

•   .   Xn 

A 

=  diag(X,, 

'^2' 

.       ,       . 

,    ^n) 

where  the  overbar  denotes  the  matrix  complex  conjugate  operator.  Based  on  this  particular 
representation  and  the  fact  that  the  state  matrices  of  Eq.  (2.2.6)  are  symmetric,  by  proper 
normalization  of  the  eigenvectors,  the  cross-orthogonality  relationship  associated  with  the 
system  are  given  by 
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(2.2.7) 
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(2.2.8) 


Equations  (2.2.7)  and  (2.2.8),  respectively,  clearly  imply  the  following  relations 
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(2.2.9) 


VA 
V 


D 
K 


K 
0 


VA 

V 


(2.2.10) 


-iT 


where  [  ]   denotes  the  complex  conjugate  transpose  operator  (  [  ]     =  )•  Note  that, 

contrary  to  Eqs.  (2.2.7-8),  no  particular  normalization  of  the  eigenvector  is  needed  in  order 
for  Eqs.  (2.2.9-10)  to  be  satisfied. 

Another  state  space  representation  of  Eq.  (2.2.3)  is  given  by 
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(2.2.11) 


Based  on  the  same  argument  discussed  earlier,  the  cross-orthogonality  conditions  that  arise 
from  this  representation  are 
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(2.2.12) 
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Again,  the  following  two  relationships  follow. 
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(2.2.13) 
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(2.2.14) 


(2.2.15) 


2.2.1     Undamped  Models 

In  the  modelling  of  structures,  it  is  often  assumed  that  the  damping  is  negligible  and 
hence  is  set  to  zero.  For  these  type  of  system  models,  the  eigenvalue  problem  in  matrix  form 
is  given  by 

MVA2  +  KV  =  [0]  (2.2.16) 

The  matrices  V  and  A  are  as  defined  earlier.  For  undamped  systems,  the  eigenvalues  X^  are 
pure  imaginary  and  the  eigenvectors  vj  are  real.  Note  that  the  eigenvalues  are  related  to  the 
system  natural  frequency  by 

Xf  =  -  CO?  (2.2.17) 

Furthermore,  by  proper  normalization  of  the  eigenvectors,  the  cross-orthogonality  relations 
associated  with  this  type  of  systems  are 
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V^MV  =  Inxn  (2.2.18) 

V'^KV  =  diag((joJ,  co^,    .  .  .  ,  oo^)  (2.2.19) 

2.2.2       Proportionally  Damped  Models 

When  the  damping  of  the  structure  is  accounted  for,  it  is  sometimes  modelled  to  be 
proportional, 

D  =  aM  +  PK  (2.2.20) 

where  a  and  P  are  real  scalars.  The  eigenvectors  of  a  proportionally  damped  system  are  the 
same  as  the  eigenvectors  associated  with  the  corresponding  undamped  system.  The 
cross-orthogonality  relationships  associated  with  proportionally  damped  models  are 

VT^MV  =  Inxn  (2.2.21) 

V^DV  =  diag(2Cia)i,  2^20)2,    .  •  •  ,  2^nO)n)  (2-2.22) 

V'^KV  =  diag(coJ,  mj,    .  .  .  ,  co^)  (2.2.23) 

where  all  variables  have  the  same  definitions  as  in  the  previous  discussion. 

A  more  detailed  development  of  the  eigenvalue  problem  can  be  found  in  the  books  of 
Meirovitch  (1986,  1980)  and  Inman  (1989). 

2.3    Experimental  Modal  Analysis 

The  vibration  characteristics  of  structures  can  also  be  measured  experimentally.  Recall 
that  the  structural  vibration  characteristics  are  given  by  natural  frequencies,  damping  ratios 
and  mode  shapes  (eigenvectors).  The  process  of  measuring  the  modal  parameters  is  known 
as  modal  testing  or  experimental  modal  analysis.  An  introductory  treatment  of  the  concept 
of  experimental  modal  analysis  can  be  found  in  the  book  of  Inman  (1989);  a  more  rigorous 
coverage  is  treated  in  the  book  of  Ewins  (1986). 

The  hardware  components  needed  in  a  modal  analysis  experiment  are  identified  in 
Figure  2.1  (Inman,  1989).  A  schematic  of  a  simple  modal  vibration  measurement  test  setup 
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is  shown  in  Figure  2.2.  Brief  descriptions  of  some  of  the  components  and  their  functions  are 
given  below. 
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T 

Ctni/^tiirf* 
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SC 

Recorder 

Signal 
Generator 

Signal  Analyzer 

Modal  Analysis 

' 

T  :=  Transducer 
SC  :=  Signal  Conditioner 


Figure  2. 1        Components  of  a  Vibration  Measurement  System 
For  Modal  Analysis  (Inman,  1989). 


The  exciter  is  used  to  deliver  the  driving  force  that  puts  the  structure  in  motion.  The 
two  most  commonly  used  exciters  are  the  shaker  and  the  impulse  hammer. 

Transducers  are  devices  that  measure  the  driving  force  as  well  as  the  response  of  the 
structure.  They  sense  mechanical  forces  or  motions  and,  then,  convert  them  into 
electrical  signals.  Mechanical  forces  are  usually  recorded  by  a  force  transducers. 
Commonly,  the  vibration  response  of  structures  is  measured  by  accelerometers. 
Accelerometers  are  mounted  directly  on  the  structure  and,  naturally,  record  the 
local  accelerations. 

Signal  conditioners  are  used  to  match  the  signal,  received  from  the  transducers,  to 
the  specification  of  the  analyzer.  These  devices  can  also  be  used  to  amplify  the 
signals. 


? 
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Figure  2.2    A  Simple  Experimental  Modal  Analysis  Setup. 

The  signal  analyzer  processes  the  electrical  signal  received  from  the  signal 
conditioners.  The  standard  type  of  analyzer  allows  time  domain  signals  to  be 
viewed  in  the  frequency  domain  via  a  Fast  Fourier  Transform  (FFT)  algorithm.  In 
an  FFT,  the  signals  are  first  filtered,  digitized  and  then  transformed  into  discrete 
frequency  spectra.  The  frequency  spectra  can  then  be  manipulated  to  compute  the 
modal  properties  of  the  structure. 

It  is  important  to  note  that  the  experimentally  identified  modal  parameters  are  usually 
affected  by  unpredictable  measurement  noise.  Typically,  natural  frequencies  are  identified 
to  within  1  to  2%  repeatability;  damping  ratios  to  within  5  to  15%  repeatability,  and  mode 
shapes  to  within  5  to  10%  repeatability.   In  practice,  due  to  testing  limitations,  the  set  of 
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structural  modal  parameters  identified  experimentally  is  incomplete  with  respect  to  the 
analytical  model.  Experimental  incompleteness  is  manifested  in  two  forms,  (i)  a  limited 
number  of  measured  modes  of  vibration,  (ii)  a  limited  number  of  eigenvector  DOF 
measurements. 

2.4    Analytical  and  Experimental  Model  Dimensions  Correlation 

One  major  prerequisite  common  to  most  model  refinement  algorithms  is  to  match  the 
number  of  degrees  of  freedom  (DOFs)  in  the  experimentally  measured  eigenvector  and  in 
the  discrete  analytical  model.  Two  different  approaches  have  been  commonly  used  to 
resolve  this  matching  problem  when  an  incomplete  set  of  eigenvector  DOFs  have  been 
measured.  The  first  approach  consists  of  reducing  the  dimension  of  the  discrete  analytical 
model  to  the  number  of  the  measured  DOFs.  The  other  approach  is  to  expand  the  measured 
eigenvector  to  the  size  of  the  analytical  model.  A  good  evaluative  survey  of  a  number  of 
analytical  model  reduction  techniques  has  been  compiled  by  McGowan  ( 1 99 1 ).  The  subject 
of  eigenvector  expansion  is  discussed  in  fair  detail  in  the  papers  by  Gysin  (1990),  and 
Zimmerman  and  Kaouk  (1992a).  In  the  next  two  sections,  some  commonly  used  model 
reduction  and  eigenvector  expansion  algorithms  are  reported  and  described. 

2.4.1     Model  Reduction  Methods 

In  this  section  the  general  framework  of  the  model  reduction  concept  is  first  presented. 
Then  formulations  of  three  commonly  used  model  reduction  techniques  are  summarized. 
Mostly,  the  concept  of  model  reduction  has  only  been  studied  for  undamped  models.  In  this 
presentation  no  attempt  has  been  made  to  generalize  these  concepts  to  damped  models. 

For  the  purpose  of  reporting  the  general  concept  of  model  reduction,  consider  the 
eigenvalue  problem  associated  with  an  undamped  model, 

MVA2  +  KV  =  [0]  (2.4.1) 
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where,  as  defined  earlier,  M  and  K  are,  respectively,  the  mass  and  stiffness  matrix;  V  is  the 
eigenvector  matrix  and  A  is  the  diagonal  eigenvalue  matrix.  Assume  that  only  a  subset  of 
the  eigenvector  DOFs  has  been  experimentally  measured. 

Equation  (2.4.1)  can  be  reordered  such  that  the  DOFs  associated  with  the  measured 
DOFs  are  in  the  upper  rows  of  the  equation. 


MoVoA^  +  KoVo  =  [0] 


(2.4.2) 


where      Vq   = 


Vu 


Ko   = 


Mo    = 


Mmm    Mmu 
Mum     Muu 


The  matrices  Mq,  Kq,  and  Vq  are,  respectively,  the  reordered  mass,  stiffness  and 
eigenvector  matrices.  In  the  above  equation,  the  subscripts  "m"  and  "u"  denote,  respectively, 
the  components  associated  with  the  measured  and  unmeasured  DOFs.  A  transformation 
matrix  P,  that  relates  matrices  Vu  and  Vm  can  be  defined  as 


Vu    =  PV, 


(2.4.3) 


A  substitution  of  this  relationship  in  the  reordered  eigenvector  matrix  Vq  results  in 


Vo   = 


V„,    =  TVr 


(2.4.4) 


Substituting  Eq.  (2.4.4)  into  Eq.  (2.4.2)  and  premultiplying  by  T^  yields  the  eigenvalue 
problem  of  the  reduced  model, 


MrVmA2  +  KrVm  =  [0]  (2.4.5) 

where      Mr  =  T'^MqT 
Kr  =  T'^KoT 

where  Mr  and  Kr  are  the  reduced  mass  and  stiffness  matrices,  respectively.  In  terms  of  the 
partitioned  matrices,  the  reduced  matrices  are  defined  as 
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Mr    =    Mmm    +    P^Mum    +    Mn,uP    +    P^MuuP 


(2.4.6) 


Kf    —    Kmni    +    P    Kum    +    KnmP    +    P    KuuP 


(2.4.7) 


At  this  point,  the  only  condition  placed  on  P  is  the  relationship  of  Eq.  (2.4.3).  Naturally, 
matrix  P  can  be  computed  directly  from  Eq.  (2.4.3)  if  the  eigenvectors  of  the  system  are 
available.  This  approach  is  known  as  the  exact  reduction  method  and  has  been  discussed  in 
the  papers  by  Kammer  (1987)  and  O'Callahan  et  al.  (1989).  The  exact  reduction  method 
requires  solving  for  a  large  number  of  eigenvectors,  which  can  be  computationally 
expensive.  The  three  reduction  methods  that  are  presented  in  the  forthcoming  discussions 
propose  alternative  techniques  to  compute  matrix  P.  The  first  two  do  not  require  the 
computation  of  the  system  eigenvalue  problem.  The  last  one  requires  the  knowledge  of  one 
eigenvalue  which  is  computationally  admissible. 

2.4.1.1    Static  Reduction 

This  reduction  method  is  often  referred  to  as  Guyan  (1965)  reduction.  In  the  static 
reduction,  the  mass  properties  associated  to  the  unmeasured  DOFs  are  assumed  negligible. 
With  that  assumption,  Eq.  (2.4.2)  can  be  written  as 


Mmm 

o' 

"Vn," 
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0 
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Vu 

v„ 


(2.4.8) 


The  second  row  of  this  matrix  equation  can  then  be  manipulated  as 


*  u    ~  ^uu  ^um  *  m 


(2.4.9) 


From  comparing  Eq.  (2.4.9)  to  Eq.  (2.4.3),  it  can  deduced  that  the  transformation 
matrix,  P,  computed  using  the  Guyan  approach  is  given  by 


—  K,„,  K 


uu   ^^um 


(2.4.10) 


The  reduced  mass  and  stiffness  matrices  can  then  be  computed  by  substituting  matrix  Pg 
for  matrix  P  in  Eqs.  (2.4.6)  and  (2.4.7).    Naturally,  the  Guyan  assumption  (Eq.  (2.4.8)) 
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suggests  that  if  the  mass  properties  of  the  omitted  DOFs  are  not  small,  the  accuracy  of  the 
Guyan  reduced  model  could  be  lacking. 

2.4.1.2    IRS  Reduction 

The  improved  reduction  method  (IRS)  was  formulated  by  O'Callahan  (1989).  It  is  an 
improvement  over  the  Guyan  reduction  in  that  it  accounts  for  the  mass  properties  of  the 
unmeasured  DOFs.  In  the  formulation  of  the  IRS  method,  the  Guyan  reduced  model  is 
corrected  to  include  the  mass  influence  of  the  unmeasured  DOFs.  This  formulation  is 
somewhat  lengthy  and  the  interested  reader  is  referred  back  to  the  paper  of  O'Callahan 
(1989)  or  the  thesis  of  McGowan  (1991).  The  transformation  matrix  P  computed  using  the 
IRS  reduction  is 


PiRS    =  Pg    +    K-u'[  Mum  +  MuuPg]  M-'Kr^ 


(2.4.11) 


The  reduced  IRS  model  is  then  computed  by  substituting  matrix  Pirs  in  Eqs.  (2.4.6-7). 

2.4.1.3    Dynamic  Reduction 

The  dynamic  reduction  was  proposed  as  another  improvement  to  the  Guyan  reduction 
(Kidder,  1973;  Miller,  1980;  and  Paz,  1984).  This  reduction  utilizes  the  dynamic  equation 
associated  with  a  single  mode  of  vibration  to  compute  the  transformation  matrix  P.  In  this 
technique,  the  transformation  matrix  P  is  arrived  at  by  considering  the  reordered  dynamic 
equation  associated  with  the  i*''  mode. 


(2.4.12) 


where  X,  is  the  i*  eigenvalue;  y^  and  Vy  are,  respectively,  the  measured  and  unmeasured 
eigenvector  DOFs  associated  with  the  i^^  mode  of  vibration.  Based  on  this  partition,  two 
equations  can  be  generated.  By  using  the  rows  that  correspond  to  the  unmeasured  DOFs 
(lower  rows),  the  following  relationship  is  obtained, 


^j  J^mm    1    ^mm 
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Yu,    =    -  [Muu?^?  +  Kuu]     '[MumXf  +  Kum]  v„^  (2.4.13) 

Hence,  the  transformation  matrix  associated  with  the  i^'^  mode  is  defined  by 


Pd.    -    - 


M„a?  +  K„„l     [m„„.X.2  +  K, 


(2.4.14) 


The  reduced  mass  and  stiffness  matrices  are  then  computed  by  using  Eqs.  (2.4.6)  and  (2.4.7). 
Note  that  different  modes  will  result  in  different  transformation  matrices  Pj  and,  hence, 

different  reduced  mass  and  stiffness  matrices. 

2.4.2     Eigenvector  Expansion  Methods 

Alternatively,  the  dimension  of  the  measured  eigenvectors  can  be  correlated  to  the 
dimension  of  the  analytical  model  by  using  eigenvector  expansion  algorithms.  The  common 
basis  of  these  algorithms  is  the  interpolation  of  the  unmeasured  eigenvector  components.  In 
the  forthcoming  sections,  two  eigenvector  expansion  algorithms  are  reviewed. 

2.4.2. 1    Dynamic  Expansion 

The  dynamic  expansion  technique  (Berman  and  Nagy,  1983)  is  one  of  the  most 
commonly  used  eigenvector  expansion  algorithms.  A  slight  modification  of  the  Berman  and 
Nagy  formulation  is  presented  here  to  accommodate  damped  systems  (Fuh  et  al.,  1984).  In 
the  formulation  of  the  dynamic  expansion,  it  is  assumed  that  the  measured  modes  satisfy  the 
eigenvalue  problem  involving  the  property  matrices  of  the  original  model  (M,  D,  K) .  For  the 
i^^  measured  mode,  this  assumption  takes  the  form 

(Xi.M  +  Xe.D  +  K)  Ye.  =  0  (2.4.15) 

where  X,e  and  Vg  are,  respectively,  the  i'*^  experimental  eigenvalue  and  eigenvector.  The 
matrices  M,  D  and  K  have  the  same  definitions  as  in  the  earlier  sections.  Assume  that  only 
a  subset  of  the  DOFs  of  eigenvector  Vg.  has  been  measured.  Equation  (2.4.15)  can  be 
reordered,  as  in  Section  2.4. 1 ,  such  that  the  measured  eigenvector  DOFs  reside  in  the  upper 
half  of  the  equation. 
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(2.4.16) 


where  Vg^and  Vgy  are,  respectively,  the  measured  and  unmeasured  DOFs  of  Vg .  The 
subscripts  "m"  and  "u"  denote  measured  and  unmeasured  components.  A  rearrangement  of 
Eq.  (2.4.16)  yields 


■^Mn 
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(2.4.17) 


From  Eq.  (2.4.17),  two  matrix  equations  in  function  of  Vg^  and  Vgy  can  be  generated. 
By  using  the  equation  associated  with  the  second  row  of  the  partition,  the  unmeasured 
components  of  the  i^^  experimental  eigenvector  are  computed  to  be 


Xeu: 


T-1 


M„„X,p..  +  D„u?le,  +  Kuu  Mum^e.  +  Dum^e.  +  K, 


'-urn 


Vem.       (2-4.18) 


Note  that  this  expansion  works  on  a  single  mode  at  a  time.  Also,  notice  that  it  involves 
the  original  analytical  model  (M,  D,  K).  This  implies  that  the  accuracy  of  the  expansion 
depends  on  the  validity  of  the  original  analytical  model. 

2.4.2.2    Orthogonal  Procrustes  Expansion 


Another  expansion  process  that  has  shown  great  promise  is  the  orthogonal  Procrustes 
(OP)  expansion  method  presented  in  the  papers  of  Smith  and  Beattie  ( 1 990)  and  Zimmerman 
and  Kaouk  (1992a).  The  technique  uses  the  general  mathematical  framework  of  the 
orthogonal  Procrustes  problem  (Golub  and  Van  Loan,  1983). 

Let  Vem  be  the  experimentally  measured  eigenvector  component  matrix  and  Vam  be  the 
corresponding  analytical  eigenvector  component  matrix.  The  essence  of  the  OP  expansion  is 
to  find  a  unitary  transformation  matrix  Pop  that  closely  rotates  Vam  into  Vgm-  This  is 
attempted  by  solving  the  following  problem. 
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Minimize 


V. 


*^amrop  Up 


(2.4.19) 


pT  p         =    T 


subject  to 

The  solution  to  this  general  problem  is  discussed  in  the  book  of  Golub  and  Van  Loan  and 
is  given  by 

Pop   =  YZT'  (2.4.20) 

where  Y  and  Z  are,  respectively  the  left  and  right  singular  matrices  of  matrix  O  defined  by 


O    =   vLVem  (2.4.21) 

Let  Vau  be  the  eigenvector  matrix  associated  with  the  unmeasured  DOFs.  In  the 
orthogonal  Procrustes  expansion,  it  is  assumed  that  the  transformation  matrix  computed 
above  also  rotates  Vau  i"to  ^^he  unmeasured  "experimental"  eigenvector  component  matrix, 
Veu,  as 

Veu    =    VauPop  (2.4.22) 

Two  different  approaches  have  been  defined  to  generate  the  expanded  experimental 
eigenvector  matrix.  The  first  is  suggested  in  the  paper  by  Smith  and  Beattie  (1990), 


V.   = 


V. 


op 


(2.4.23) 


In  this  definition,  the  expanded  "experimental"  eigenvector  matrix  is  the  rotated  analytical 
eigenvector  matrix.   The  actual  experimental  measured  eigenvector  components  Vgm  are 
replaced  by  matrix  VamPop.     The  advantage  of  the  approach  is  that  the  resulting 
"experimental"  eigenvector  satisfies  the  cross-orthogonality  conditions  (see  Section  2.4). 
The  second  viewpoint  was  proposed  by  Zimmerman  and  Kaouk  (1992a), 

V. 


v.  = 


em 


V    P 

*  au*^  op 


(2.4.24) 


Here,  the  unaltered  eigenvector  components  measured  experimentally  are  inserted  in  Vg.  In 
this  viewpoint,  if  cross-orthogonalization  of  the  expanded  experimental  eigenvectors  is 
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required,  a  separated  orthogonality  algorithm  can  be  used  (Section  2.5).  In  the  paper  by 
Smith  at  al.  (1993),  it  is  shown  that  for  actual  model  refinement  problems,  both  viewpoints 
give  equivalent  eigenvector  expansion  results.  However,  for  damage  detection  problems, 
a  preliminary  study  indicates  that  eigenvectors  expanded  by  using  the  viewpoint  defined  in 
Eq.  (2.4.24)  give  better  assessment  of  the  damage. 

2.5    Eigenvector  Orthogonalization 

Most  matrix  update  algorithms  require  the  measured  eigenvectors  to  satisfy  a 
cross-orthogonality  condition.  This  is  especially  true  in  matrix  update  algorithms  in  which 
(i)  the  model  of  the  structure  is  assumed  undamped  and  the  modelling  errors  are  assumed  to 
be  in  only  one  of  the  two  property  matrices  (M  or  K  is  assumed  correct)  or  (ii)  the  system  is 
modeled  by  a  proportionally  damped  model  with  errors  in  only  two  of  the  three  property 
matrices  (M,  D  or  K).  In  these  situations,  in  order  to  insure  symmetry  of  the  updated  property 
matrices,  it  is  required  that  the  expanded  experimental  eigenvectors  be  orthogonal  with 
respect  to  the  property  matrix  assumed  accurate.  This  situation  is  encountered  in  the  optimal 
update  algorithms  developed  by  Baruch  and  Bar  Itzhack  (1978),  Kabe  (1985),  Kammer 
(1985),  and  Smith  and  Beattie  (1991),  Zimmerman  and  Kaouk  (1992b),  Kaouk  and 
Zimmerman  (1993b)  among  others.  In  most  of  these  algorithms,  it  is  assumed  that  the  mass 
matrix  is  correct.  This  assumption  is  used  in  a  number  of  the  model  refinement  algorithms 
since  the  inertial  properties  of  structures  are  known  to  a  good  extent.  In  these  cases,  one 
would  expect  the  expanded  experimental  eigenvectors  to  be  mass  orthogonal.  However, 
because  of  measurements  errors,  this  condition  rarely  occurs.  For  this  reason,  a  great  deal  of 
effort  was  focussed  on  the  development  of  mass  orthogonalization  techniques.  Some  of  the 
most  notable  work  in  that  area  was  performed  by  Targoff  (1976),  Baruch  and  Bar  Itzhack 
(1978),  and  Baruch  (1979).  In  the  next  two  sections,  two  orthogonalization  techniques 
(Baruch  and  Bar  Itzhack,  1978;  Baruch,  1979)  are  discussed.  Both  techniques  are  mass 
orthogonalization  techniques;  however,  with  obvious  modifications,  these  techniques  can  be 
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adopted  to  solve  the  orthogonalization  problem  of  the  eigenvectors  with  respect  to  the 
stiffness  or  the  damping  matrices. 

2.5.1  Optimal  Weighted  OrthogonaUzation 

The  essence  of  the  standard  mass  orthogonalization  technique  is  to  modify  the  measured 
eigenvectors  such  that  the  mass  cross-orthogonality  condition  is  satisfied.  Baruch  and  Bar 
Itzhack  ( 1 978)  proposed  an  elegant  solution  to  that  problem.  An  overview  of  their  problem 
statement  and  solution  is  given  below. 

Assuming  that  Vg  is  a  matrix  of  expanded  experimental  eigenvectors  that  need  to  be 
mass  orthogonalized.  The  present  formulation  searches  for  the  optimal  mass  weighted 
change  of  matrix  Vg  such  that  the  mass  cross-orthogonality  condition  is  satisfied.  This 
problem  is  cast  as 

Minimize      ||  N  (Veo  -  Vg)  ||f  (2.5.1) 

subject  to      Vjo  M  Veo  =  I  (2-5.2) 

where  N=M'^2  and  M  is  the  mass  matrix.  By  means  of  a  Lagrange  multiplier,  Eq.  (2.5.2) 
can  be  incorporated  into  Eq.  (2.5.1);  then  the  application  of  the  optimality  conditions  yields 
the  following  expression  for  Vgo, 

Veo  =  Ve(vTMVe)''^'  (2-5.3) 

Before  being  incorporated  into  the  orthogonalization  process,  the  measured 
eigenvectors  have  to  be  unit  mass  normalized,  i.e., 

v.,  =  v,(vTMvJ""'  (2.5.4) 

where  Vg  is  the  i''^  expanded  experimental  eigenvector  (i'"^  column  of  Vg). 

2.5.2  Selective  Optimal  Orthogonalization 

Some  structures  exhibit  rigid  body  modes  (modes  with  zero  eigenvalues).  It  is  desirable 
to  preserve  these  rigid  body  modes  in  the  refinement  process.  However,  some  matrix  update 
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algorithms  require  the  rigid  body  modes  and  the  experimental  eigenvectors  to  be  mass 
orthogonal  to  insure  symmetry  of  the  updated  property  matrices  (see  Chapter  6).  Naturally, 
the  rigid  body  modes  will  be  corrupted  if  they  are  incorporated  along  with  the  expanded 
experimental  eigenvectors  in  the  above  orthogonalization  process.  Thus,  Baruch  (1979) 
presented  a  modification  of  the  procedure  in  Section  2.5. 1  to  deal  with  such  a  problem.  The 
resulting  problem  is  a  selective  orthogonalization  and  is  formulated  as  follows, 

minimize      ||  N  (Vgo  —  Vg)  ||p  (2.5.5) 

Subject  to      Vjo  M  Veo  =  I  (2.5.6) 

and      Vjo  M  Vr  =  [0]  (2.5.7) 

In  the  above  equations,  Vg  and  Vr  are,  respectively,  the  expanded  experimental  eigenvector 
matrix  and  the  rigid  body  mode  matrix.  Again,  the  Lagrange  multiplier  is  used,  and  the 
orthogonalized  experimental  eigenvector  matrix  that  satisfies  the  conditions  in  Eqs.  (2.5.6) 
and  (2.5.7)  is  found  to  be 

(2.5.8) 
where     Q  =  Vg  -  VrV?"MVe 

Note  that,  as  in  previous  process,  the  expanded  experimental  eigenvectors  have  to  be 
unit  mass  normalized. 

2.6    Load  Path  Preservation 

Many  matrix  update  algorithms  introduce  additional  load  paths  in  their  updated  models, 
i.e.,  elements  of  the  mass,  damping  or  stiffness  matrices  that  were  originally  zero  may 
become  nonzero.  Whether  or  not  preserving  the  original  load  path  is  a  practical  problem  is 
still  a  matter  of  current  debate.  It  seems  that  for  damage  assessment  of  truss  structures  it  is 
desired  to  maintain  load  paths.  In  the  paper  of  Zimmerman  and  Kaouk  (1992a),  an  iterative 
approach  to  preserve  the  load  path  of  the  original  property  matrices  was  developed.  The 


29 


approach  was  presented  in  the  context  of  the  symmetric  eigenstructure  assignment  model 
refinement  algorithm  (discussed  in  Chapter  4);  however,  its  application  can  also  be  extended 
to  other  model  refinement  algorithms.  In  Figure  2.3,  a  flow  chart  of  the  iterative  load 
preservation  algorithm  is  presented.  The  procedure  is  illustrated  for  a  general  model 
refinement  scenario  in  which  all  three  property  matrices  (M,  D,  K)  are  being  updated. 
However,  it  can  be  easily  modified  to  accommodate  other  refinement  problems. 
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Figure  2.3    Flow  Chart  of  the  Iterative  Load  Path  Preservation  Algorithm 


In  the  flow  chart,  the  matrices  Ma,m.  Da,m>  and  Ka,ni  are  respectively  the  adjusted 
masked  mass,  damping  and  stiffness  matrices  defined  by 
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Ma,m  =  Ma  O  Mm 

Da,m  =  Da  O  Dm  (2.6.2) 

Ka,m  =  Ka  O  Km 

where  Ma,  Da,  and  Ka  are  the  adjusted  mass  damping  and  stiffness  matrices.  The  matrices 
Mm,  Dm,  and  Km  are  the  masking  matrices  associated  with  the  original  mass,  damping  and 
stiffness  matrix.  By  definition,  the  masking  matrix.  Am,  associated  with  matrix  A  is  given 
by 

Am(i,j)  =  1    if  A(i,j)  ^  0 

(2.6.1) 
Am(i,j)  =  0    if  A(i,j)  =  0 

In  Eqs.  (2.6.2),  the  operator  O  is  the  element-by-element  (scalar)  matrix  multiplication.  Let 
B  and  C  be  two  nxn  matrices,  then  the  element-by-element  multiplication  of  B  and  C  is  given 
by 

S  =  BOC    ^    S(i,j)  =  B(i,j)  *  C(i,j)         i,j  =  l,...,n  (2.6.3) 

At  every  iteration,  the  norms  of  the  matrix  differences  between  corresponding  adjusted 
and  adjusted  masked  property  matrices  are  computed.  At  a  given  iteration,  if  the  three 
computed  norms  are  equal  to  zero  or  within  user  set  limits,  then  the  load  paths  of  the  original 
three  property  matrices  have  been  exactly  achieved  or  achieved  within  user  state  guidelines. 
Thus,  the  procedure  is  halted,  and  the  refined  model  consists  of  the  adjusted  property 
matrices  computed  at  that  particular  iteration.  It  should  be  noted  that  there  is  no  formal 
guarantee  of  convergence  in  using  this  iterative  procedure.  Experience  gained  in  using  the 
present  algorithm  indicates  that  if  the  experimental  modal  data  are  consistent  with  the 
sparsity  pattern,  the  procedure  will  converge.  Consistent  data  means  that  there  exist  mass, 
damping  and  stiffness  matrices  that  have  the  same  sparsity  pattern  as  the  original  matrices 
and  also  exhibit  the  measured  test  data.  Otherwise,  if  the  data  are  inconsistent,  the  original 
sparsity  patterns  will  not  be  exactly  preserved.  In  this  case,  the  added  load  path  terms  of  Ma, 
Da  and  Ka  which  should  be  zero  will  be  closer  to  zero  after  application  of  the  algorithm. 


CHAPTER  3 

INVERSE  /  HYBRID  APPROACH  FOR  FINITE 

ELEMENT  MODEL  REFINEMENT 


3.1    Introduction 

The  inverse  eigenvalue  problem  is  concerned  with  the  construction  of  the  property 
matrices  (mass,  damping  or  stiffness)  of  a  dynamic  model  using  experimentally  measured 
modal  data.  These  techniques  require  complete  modal  properties.  Thus,  for  an  n-DOF 
model,  n  natural  frequencies,  damping  ratios  and  mode  shapes  (eigenvectors)  must  be 
measured,  and  the  identified  mode  shapes  must  be  of  dimension  n.  Due  to  practical  testing 
limitations,  this  is  rarely  accomplished  for  typical  large  structures.  In  this  chapter,  the 
application  of  the  inverse  problem  is  extended  to  model  updating  by  combining 
experimental  measurements  and  original  analytical  FEM  modal  information.  Again, 
refinement  implies  correlating  the  measured  and  analytical  modal  properties. 

3.2    Theoretical  Formulation 

The  dynamic  structure  under  consideration  is  assumed  to  be  successfully  modelled  by 
an  n-DOF  proportionally  damped  nongyroscopic  and  noncirculatory  (symmetric  property 
matrices)  FEM.  The  free  vibration  motion  of  such  a  dynamic  structure  can  be  analytically 
represented  by  a  differential  equation  of  the  form 

Mw(t)  +  Dw(t)  +  Kw(t)  =  0  (3.2.1) 

where  the  variables  M,  D,  and  K  are  nxn  real  symmetric  matrix  models  of  the  mass,  damping 
and  stiffness  properties  of  the  structure.  The  nxl  time  varying  vector  w(t)  represents  the  n 
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displacements  of  the  n-DOF  model  of  the  system.  The  overdots  represent  differentiation 
with  respect  to  time.  The  eigenvalue  problem  associated  with  the  differential  equation 
shown  in  Eq.  (3.2.1)  is  of  the  form 

MVj  If  +  DVj  \  +  KVj  =  0  (3.2.2) 

where  Ai  and  Vj  are,  respectively,  the  eigenvalue  and  eigenvector  of  the  i'^  mode  of  vibration. 
In  this  problem,  it  is  assumed  that  the  accuracy  of  the  original  FEM  is  lacking  and,  hence, 
needs  improvement.  Furthermore,  it  is  assumed  that  the  inaccuracy  of  the  original  FEM  is     . 
solely  due  to  modeling  errors  in  the  stiffness  and  damping  properties. 

The  model  refinement,  proposed  herein,  exploits  the  cross-orthogonality  relations  that 
arise  from  the  symmetric  nature  of  the  property  matrices  and  the  proportional  damping 
assumption.  As  discussed  in  Chapter  2,  by  proper  normalization  of  the  eigenvectors  these 
cross-orthogonality  relations  have  the  form 

VT  M  V  =  Inxn  (3.2.3a) 

VT  D  V  =  diag(2Ci(0i,    .  .  .  ,2Cn(0n)  =  2  (3.2.3b) 

V'^KV  =  diag(a)i2,    .  .  .  .con^)  =  Q  (3.2.3c) 


V  =  [    V,    ,    ...    ,    Vn 

where  tOj  and  ^j  are  the  natural  frequency  and  damping  ratio,  respectively,  of  the  i"^  mode 
of  the  structure.  The  matrix  Inxn  is  the  nxn  identity  matrix.  It  is  important  to  recognize  that 
Eq.  (3.2.3a)  represents  necessary  and  sufficient  conditions  for  conserving  symmetry  and 
damping  proportionality  when  updating  the  stiffness  and  damping  properties  of 
proportionally  damped  systems. 

Suppose  that  p  (p  «  n)  modes  of  an  existing  structure  have  been  experimentally 
identified  (mode  shapes  or  eigenvectors,  frequencies  and  damping  ratios).  Assume  that  the 
dimension  of  measured  eigenvectors  is  equal  to  the  dimension  of  the  FEM,  i.e.,  all  n 
components  of  the  measured  eigenvectors  are  available.  It  is  widely  accepted  that  in  the 


absence  of  specific  experimental  measurements  a  good  approximation  to  the  unmeasured 
modes  is  their  corresponding  analytical  modal  information.  With  that  in  mind,  a  complete 
hybrid  set  of  modal  data  is  generated  by  combining  experimental  and  analytical  information 

'ea~l''e        »aJ 


(3.2.4) 


f^ 

"Qe    0 

^ea  = 

0      ^a 

2ea  = 

"Se    0    " 
0     2a 

where  V  is  the  eigenvector  matrix;  Q.  and  Z  are  diagonal  matrices  of  frequencies  squared  and 
damping  ratios,  respectively.  The  subscripts  e  and  a  denote,  respectively,  experimental  and 
analytical  sets.  At  this  point,  the  complete  "hybrid"  set  of  eigenvectors,  Vea,  does  not  satisfy 
the  cross-orthogonality  conditions  defined  in  Eq.  (3.2.3a);  thus  the  conditions  in  Eqs. 
(3.2.3b)  and  (3.2.3c)  are  not  met.  One  possible  solution  to  this  problem  is  to  modify  all  of 
Vea  in  an  optimal  way  to  comply  with  the  orthogonality  requirement.  This  approach  treats 
all  parts  of  Vga  equally,  and  thus  overlooks  the  fact  that  the  experimental  modes  are  known 
with  a  higher  confidence.  Naturally,  a  technique  that  assigns  a  higher  credibility  to  the 
experimental  eigenvectors,  Vg,  is  preferable.  This  can  be  achieved  by  incorporating  Vga  into 
the  orthogonalization  process  group  by  group  in  the  order  of  their  descending  credibility 
(experimental  then  analytical).  If  the  problem  is  set  such  that  the  experimental  modes,  Ve, 
are  corrected  first,  it  is  clear  that  the  analytical  modes,  Va,  will  be  subject  to  larger  correction 
when  incorporated  into  the  orthogonalization  process  since  they  will  be  subject  to  more 
constraints. 

The  experimental  modes,  Vg,  are  orthogonalized  by  using  the  orthogonalization 
technique  formulated  by  Baruch  and  Bar  Itzhack  (1978).  The  general  formulation  of  this 
orthogonalization  technique  is  discussed  in  Chapter  2.  For  this  particular  application,  the 
problem  consists  of  finding  the  matrix  ,Veo,  that  satisfies  the  cross-orthogonality  condition, 
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VL  M  Veo  =  Ipxp  (3.2.5) 

and  that  minimizes  the  weighted  Euclidean  norm, 


(3.2.6) 


^  =||N(Veo-  Ve)||F 

where        N  =  M''^ 
The  solution  to  this  problem,  as  reported  in  Chapter  2,  is 

Veo  =  Ve(vlMVe)"'^'  (3-2.7) 

The  next  step  is  to  invoke  the  orthogonality  requirement  on  the  analytical  eigenvector 
matrix,  Va,  by  searching  for  a  matrix  Vao  that  satisfies  the  following  two  conditions, 

Vio  M  Vao  =  I(n-p)x(n-p)  (3.2.8) 

and      Vio  M  Veo  =  [0]  (3.2.9) 

while  minimizing  the  objective  function, 

r=||N(Vao-  Va)||F 

(3.2.10) 
where        N  =  M'/'^ 

A  similar  problem  was  also  treated  by  Baruch  (1980)  in  a  different  context.  A  brief 
discussion  of  the  solution  approach  can  be  found  in  Chapter  2.  The  set  of  eigenvectors,  Vao 
that  satisfies  this  problem  is  given  by 

Vao  =  Q(Q''MQ)"'^^ 

(3.2.11) 
Q  =  Va  -  VeoVjoMVa 

Clearly,  the  resultant  matrix,  Vgao  =  [  Vgo  ,  Vao  ]>  satisfies  Eq.  (3.2.3a).  The  corrected 
stiffness  and  damping  matrices  are  then  computed  using  Eqs.  (3.2.2b)  and  (3.2.2c). 
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Ka  =  MVeao  Q  eaVjaoM  (3.2.12a) 

Da  =  MVeao  ©  eaV^M  (3.2.12a) 

where       ©ea  =  2neaQea 

Note  that  the  matrices  Ka  and  Da  computed  from  Eqs.  (3.2.12a,b)  will  be  symmetric. 
The  above  formulation  suggests  that  the  system  modelled  by  the  original  mass  matrix  (M) 
and  the  stiffness  (Kg)  and  damping  (Dg)  matrices  computed  from  Eqs.  (3.2.12a,b)  will  have 
eigenvectors  Veao,  frequencies  Dea.  and  damping  ratios  riea-  Some  structures  exhibit  rigid 
body  modes  of  vibration.  Commonly,  it  is  desirable  to  preserve  these  rigid  body  modes  in  the 
updated  model.  The  above  formulation  also  suggest  that  the  updated  model  will  preserve  the 
original  rigid  body  modes. 

The  procedure  developed  above  can  be  easily  contracted  to  address  the  case  when  the 
system  model  does  not  account  for  the  effects  of  damping  (undamped  model).  The 
contraction  can  be  obtained  by  setting  to  zeros  matrices  D,  ©gg,  and  flga  in  Eqs.  (3.2.2^). 

The  computational  burden  limits  the  size  of  the  FEM  which  can  be  updated  (order  of 
200).  Essentially,  the  limiting  factor  is  that  all  mode  shapes  of  the  structure  that  are  not 
available  from  experimental  measurements  must  be  calculated  analytically  by  solving  the 
eigenvalue  /  eigenvector  problem. 

3.3    Numerical  Illustration 

The  system  addressed  in  this  investigation  is  the  commonly  used  eight  degrees  of 
freedom  model  shown  in  Figure  3.1.  This  model  was  developed  by  Kabe  (1985)  to  give  a 
common  testbed  for  the  evaluation  of  the  performance  of  model  refinement  algorithms.  An 
original  undamped  analytical  model  of  the  system  was  generated  by  using  the  mass  and 
stiffness  properties  shown  in  Figure  3.1.  The  elements  of  the  original  analytical  stiffness 
matrix  are  displayed  in  the  second  column  of  Table  3. 1 .  In  this  problem,  it  is  assumed  that  the 
original  stiffness  matrix  of  the  model  is  incorrect.  These  inaccuracies  were  simulated  by 
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mi  =  0.001      m8  =  0.002      mj  =  1.0     j  =  2, ... ,  7 
ki  =  1000     k2=10     k3  =  900     k4=100     kg  =1.5     k6  =  2.0 

Figure  3.1    Kabe's  Problem:  Analytical  Test  Structure. 


using  incorrect  stiffness  constants  for  most  of  the  springs.  The  elements  of  the  actual  correct 
stiffness  matrix  are  given  in  the  last  column  of  Table  3.1.  Note  that  the  present  model 
refinement  problem  is  very  challenging  and  because  of  the  large  difference  between  the 
stiffness  matrices  of  the  original  and  exact  model.  In  this  study,  two  cases  will  be  considered. 
In  the  first  case,  it  is  assumed  that  only  the  modal  parameters  (eigenvalue  and  eigenvector)  of 
the  first  mode  were  measured.  In  the  other,  the  modal  parameters  of  first  three  modes  are 
assumed  to  be  available.  In  both  cases,  the  measured  eigenvectors  are  supposed  to  be  full  (all 
degrees  of  freedom  of  the  eigenvector(s)  are  measured). 

The  main  objective  of  this  investigation  is  to  compare  the  performance  of  the 
inverse/hybrid  method  to  the  algorithm  proposed  by  Baruch  and  Bar  Itzhack  (1978).  The 
Baruch  and  Bar  Itzhack  model  update  technique  is  one  of  the  most  commonly  used  model 
refinement  algorithm.  One  of  the  main  reasons  for  its  common  use  is  because  the  algorithm 
is  computationally  efficient.  The  updated  stiffness  matrices  generated  by  using  the  Baruch 
and  Bar  Itzhack  approach  for  the  one  mode  and  three  mode  cases  are  shown  in  the  third  and 
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fourth  column  of  Table  3.1,  respectively.  The  fifth  and  sixth  columns  of  Table  3. 1  display  the 
results  of  using  the  inverse/hybrid  algorithm  for  the  one  mode  and  three  modes  cases.  For 
both  cases,  it  is  clear  that  the  performances  of  both  algorithms  is  lacking  in  predicting  the 
exact  stiffness  matrix.  As  expected.  It  can  also  be  seen  that  both  algorithms  generate  better 
results  when  using  three  measured  modes.  A  comparison  of  the  results  generated  using  the 
Baruch  and  Bar  Itzhack  approach  and  the  inverse/hybrid  method  shows  that  both  algorithms 
give  the  same  type  of  results.  This  similarity  of  results  was  also  encountered  in  other 
example  problems. 

3.4    Summary 

A  model  refinement  approach  that  uses  a  hybrid  set  of  experimental  and  analytical 
modal  properties  was  formulated.  The  developed  approach,  termed  the  inverse/hybrid 
algorithm,  was  illustrated  by  using  a  computer  simulated  example.  Part  of  the  evaluation  of 
the  proposed  algorithm  was  the  comparison  of  its  performance  with  the  performance  of  the 
Baruch  and  Bar  Itzhack  algorithm.  It  was  found  that  both  algorithms  give  the  same  type  of 
results.  However,  the  computations  involved  in  the  inverse/hybrid  approach  exceed  those 
involved  in  the  Baruch  and  Bar  Itzhack  algorithm.  Essentially,  the  inverse/hybrid  approach 
requires  the  computation  of  all  eigenvalues  and  eigenvectors  of  the  structure  that  are  not 
available  from  experimental  measurements.  As  will  be  illustrated  in  the  forthcoming 
chapter,  the  Baruch  type  approaches  are  not  suited  for  damage  assessment  applications.  For 
these  reasons,  further  development  of  the  Inverse/Hybrid  algorithm  was  not  investigated  and 
new  formulations  (Chapters  4  &  6)  were  considered. 
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Table  3. 1    Kabe's  Problem:  Elemental  Stiffness  Components. 


Baruch 

Inverse/Hybrid 

Element  # 

Original 

Mode  1 

Modes  1-3 

Mode  1 

Modes  1-3 

Exact 

(1,1) 

2.0 

2.0 

2.0 

2.0 

2.0 

1.5 

(1,2) 

-2.0 

-2.0 

-3.0 

-2.0 

-3.0 

-1.5 

(1,3) 

0.0 

0.0 

-0.1 

0.0 

-0.1 

0.0 

(1,4) 

0.0 

0.0 

-0.1 

0.0 

-0.1 

0.0 

(1,5) 

0.0 

0.0 

-0.1 

0.0 

0.1 

0.0 

(1,6) 

0.0 

0.0 

0.1 

0.0 

0.1 

0.0 

(1,7) 

0.0 

0.0 

0.1 

0.0 

0.0 

0.0 

(1,8) 

0.0 

0.0 

0.0 

0.0 

0.0 

0.0 

(2,2) 

1512.0 

1508.6 

1024.2 

1510.6 

1024.3 

1011.5 

(2,3) 

-10.0 

-31.5 

-68.5 

-21.5 

-71.4 

-10.0 

(2,4) 

0.0 

-8.9 

-9.0 

-15.0 

-8.4 

0.0 

(2,5) 

0.0 

-8.9 

-20.9 

-15.0 

23.8 

0.0 

(2,6) 

0.0 

-21.6 

38.5 

-11.5 

35.2 

0.0 

(2,7) 

0.0 

-3.9 

-9.1 

-1.6 

-8.7 

0.0 

(2,8) 

0.0 

-0.1 

0.1 

-0.1 

0.1 

0.0 

(3,3) 

1710.0 

1574.6 

1560.8 

1624.1 

1612.2 

1110.0 

(3,4) 

0.0 

^4.9 

^9.8 

-75.6 

-76.2 

0.0 

(3,5) 

-200.0 

-244.8 

-244.1 

-275.6 

-276.9 

-100.0 

(3,6) 

0.0 

-136.5 

-123.4 

-86.2 

-74.5 

0.0 

(3,7) 

0.0 

-24.6 

50.8 

-13.1 

46.8 

0.0 

(3,8) 

0.0 

-0.4 

-0.1 

-0.4 

0.0 

0.0 

(4,4) 

850.0 

1083.1 

1087.9 

1102.1 

1099.9 

1100.0 

(4,5) 

-200.0 

3.27 

25.1 

51.8 

48.0 

-100.0 

(4,6) 

-200.0 

-254.3 

-242.2 

-276.5 

-274.9 

-100.0 

(4,7) 

0.0 

-10.1 

27.0 

-17.2 

30.8 

0.0 

(4,8) 

0.0 

-0.2 

-0.1 

-0.2 

-0.2 

0.0 

(5,5) 

850.0 

1082.2 

1089.5 

1101.4 

1101.1 

1100.0 

(5,6) 

0.0 

^5.2 

-49.5 

-76.5 

-75.1 

0.0 

(5,7) 

0.0 

-10.1 

-12.4 

-17.2 

-11.7 

0.0 

(5,8) 

0.0 

-0.2 

-0.1 

-0.2 

-0.2 

0.0 

(6,6) 

1714.0 

1576.4 

1565.0 

1617.5 

1610.8 

1112.0 

(6,7) 

-10.0 

-34.8 

-82.6 

-23.1 

-76.2 

-10.0 

(6,8) 

^.0 

^.4 

-A.l 

-AA 

-4.5 

-2.0 

(7,7) 

1512.0 

1507.5 

1027.9 

1510.1 

1028.0 

1011.5 

(7,8) 

-2.0 

2.1 

-4.1 

-2.1 

^.2 

-1.5 

(8,8) 

6.0 

6.0 

6.0 

6.0 

6.0 

3.5 

\..F 


CHAPTER  4 

SYMMETRIC  EIGENSTRUCTURE  ASSIGNMENT  MODEL 

REFINEMENT  ALGORITHM 


4.1    Introduction 

Eigenstructure  assignment  is  a  control  concept  used  to  alter  the  transient  response  of 
linear  systems.  This  is  done  by  forcing  the  system  to  have  some  predetermined  eigenvalues 
and  eigenvectors.  A  detailed  overview  of  eigenstructure  assignment  theories  can  be  found  in 
the  paper  by  Andry  et  al.  (1983).  Inman  and  Minas  (1990),  Zimmerman  and  Widengren 
(1989, 1990),  and  Widengren  (1989)  have  developed  model  refinement  algorithms  based  on 
the  mathematical  framework  of  eigenstructure  assignment.  The  basic  idea  of  these  model 
refinement  techniques  is  to  design  the  pseudo-controller  which  is  required  to  produce  the 
measured  modal  properties  (natural  frequencies,  damping  ratios  and  mode  shapes)  with  the 
original  finite  element  model  (FEM)  of  the  structure.  The  pseudo-controller  is  then 
translated  into  matrix  adjustments  applied  to  the  initial  FEM. 

In  this  work,  the  eigenstructure  assignment  based  model  refinement  algorithm  proposed 
by  Zimmerman  and  Widengren  (1989,  1990)  is  extended  to  better  approach  the  damage 
assessment  problem.  A  subspace  rotation  algorithm  is  developed  to  enhance  eigenvector 
assignability.  Finally,  the  enhanced  algorithm  is  tested  and  compared  to  other  techniques  on 
both  "simulated"  and  actual  experimental  data. 

4.2    Problem  Formulation 

In  this  section,  a  review  of  the  Zimmerman  and  Widengren  (1989,  1990)  refinement 
technique,  which  is  termed  the  symmetric  eigenstructure  assignment  model  refinement 
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algorithm  (SEAMRA),  is  presented.  This  review  is  essential  in  order  to  properly  introduce 
and  discuss  the  extension  and  improvement  proposed  in  this  work. 

4.2.1     Standard  Eigenstructure  Assignment  Formulation 

Consider  the  standard  differential  equation  of  motion  of  an  n  degrees  of  freedom 
damped,  nongyroscopic  and  noncirculatory  structure  with  control  feedback, 

Mw(t)  +  Dw(t)  +  Kw(t)  =  Bou(t)  (4.2.1) 

Again,  M,  D,  and  K  are  n-by-n  real  symmetric  matrix  models  of  the  mass,  damping  and 
stiffness  properties  of  the  structure.  Assume  that  these  matrices  were  generated  using  the 
finite  element  method.  The  nxl  time  varying  vector,  w(t),  represents  the  n  displacements 
of  the  n-DOF  FEM  of  the  system.  The  overdots  represent  differentiation  with  respect  to 
time.  In  control  terminology,  Bq  is  the  nxm  (m  «  n)  control  influence  matrix  describing 
the  actuator  force  distributions  and  u(t)  is  the  mx  1  vector  of  output  feedback  control  forces 
defined  by 

u(t)  =  Fy(t)  (4.2.2) 

In  Eq.  (4.2.2),  F  is  the  mxr  feedback  gain  matrix  and  y(t)  is  the  rxl  output  of  sensor 
measurements  defined  by 

y(t)  =  CqW  +  CjW  (4.2.3) 

in  which  Cq  and  Ci  are  the  rxn  output  influence  matrices  corresponding  to  position  and 
velocity,  respectively.  A  substitution  of  Eqs.  (4.2.2)  and  (4.2.3)  into  Eq.  (4.2. 1)  yields 

Mw(t)  +  (D  -  BoFC,)w(t)  +  (K  -  BoFCo)w(t)  =  0  (4.2.4) 

It  is  clear,  from  Eq.  (4.2.4),  that  the  feedback  controller  results  in  residual  changes,  BqFCo 
and  BqFCi,  to  the  stiffness  and  damping  matrices,  respectively.  These  changes  can  be 
viewed  as  perturbations  to  the  initial  finite  element  model  (FEM)  such  that  the  adjusted  FEM 
matches  closely  the  experimentally  measured  modal  properties.  The  adjusted  FEM  consists 
of  the  original  mass  matrix  and  the  adjusted  stiffness  and  damping  matrices  given  by 
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Kg     —     K     —     BqFCq 


Da  =  D  -  BoFC, 


(4.2.5) 


Assume  that  modal  analysis  of  the  structure  under  consideration  has  been  performed  and  that 
p  modes  (p  eigenvalues  Xg.,  and  p  eigenvectors  Vg )  have  been  identified.  As  discussed  earlier 
in  Chapter  2,  in  practice  p  is  typically  much  less  than  n.  The  feedback  gain  matrix  F,  such 
that  the  adjusted  FEM  eigendata  matches  the  experimental  modal  parameters,  is  computed 
using  standard  eigenstructure  assignment  theories  (Andry  et  al.  1983): 


F  =  [Z  -  AuV]<^ 


[Co    C,] 


w    w 

WA  WA 


-.  >.  -1 


(4.2.6) 


where 


A  = 


0 


I 


nxn 


-M-^K     -M-'D 


B  = 


0 

m-'k 


T  =  [B    P] 


A  =  T-'AT  = 


■Au' 


B  =  T-'B  = 


0 


V  =  T" 


w     w 

WA    WA 


A  =  diag(Xe,,  Xe,, ...,  XeJ  W  =  [v,^^,y^^, ...,  VeaJ 


z  =  s 


w 

WA 


w 

WA 


1  _ 


The  overbar  in  the  above  equations  indicates  the  complex  conjugate  operator.  The  vectors 
Vga.  in  matrix  W  are  the  expanded  "best  achievable"  eigenvectors  associated  with  the 
experimentally  measured  eigenvectors  Vg,     An  explanation  of  the  concept  of  "best 
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achievable"  eigenvectors  is  discussed  in  Section  4.2.3.   The  submatrix  P  of  matrix  T  is 
arbitrary  as  long  as  T  is  invertible. 

At  this  point,  the  variables  Bq,  Cq  and  Cj  are  still  arbitrary.  A  random  selection  of  these 
variables  will  usually  result  in  nonsymmetric  perturbation  matrices  and,  consequently 
nonsymmetric  adjusted  stiffness  and  damping  matrices.  This  clearly  conflicts  with  the 
fundamental  symmetry  requirement  of  most  structures'  FEM.  In  the  formulation  of  Inman 
and  Minas  (1990),  the  resulting  perturbation  matrices  from  the  pseudo-controller  are  forced 
to  be  symmetric  through  a  nonlinear  unconstrained  optimization  problem.  Zimmerman  and 
Widengren  (1989,  1990)  proposed  a  non-iterative  and  computationally  more  efficient 
approach  to  satisfy  the  symmetry  requirement.  This  approach  in  discussed  in  the  following 
section. 

4.2.2     Symmetric  Eigenstructure  Assignment  Formulation 

The  perturbation  matrices  are  symmetric  if  the  following  conditions  are  met, 

(4.2.7) 


BoFCo  =  CjE^Bj  (a) 


BqFCi  =  C|fTbJ  (b) 

At  this  point,  two  additional  assumptions  are  made.  As  a  prerequisite  to  the  existence  of  the 
inverse  of  some  matrices  used  in  the  computations,  it  is  assumed  that  the  number  of  pseudo 
sensors  and  actuators  is  equal  to  twice  the  number  of  measured  modes  (m=r=2p).  The  other 
assumption  consists  of  restricting  the  matrices  Cq  and  Ci  by  the  conditions 

Co  =  GoBj  (a) 


(4.2.8) 

0 


Ci  =  G,Bj  (b) 


where  Gq  and  Gi  are  mxm  invertible  matrices.  A  substitution  of  Eq.  (4.2.8)  into  Eq.  (4.2.7) 
simplifies  the  symmetry  conditions  to  the  following  relationships. 


FGq  =  GjF'f"  (a) 

FGi  =  gIfT  (b) 


(4.2.9) 
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By  using  the  conditions  in  Eq.  (4.2.9),  along  with  the  expression  for  the  feedback  gain  matrix 
(Eq.  (4.2.6)),  a  necessary  but  not  sufficient  condition  on  Gq  and  G],  for  symmetric 
perturbation  matrices,  is  expressed  in  the  form  of  a  generalized  algebraic  Riccati  equation. 


A,X  +  XA2  +  XA3X  +  A4  =  [0] 


(4.2.10) 


where 


X  =  Gf  iGo 


A,  = 


_*— 1— 1/**  \        * 

a  a    'o    Mat    -  xal  —  x 


a 


A  *  —  1  —  1         * 

A,  =  T  a    'o    'a 


A.  = 


*— 1— 1/**  \— * 

3  =  T  a    'a    Max    —  xala 


A  *  —   1  —   1  *  T 

A.  =  a  a    'a    'a    -  L 


X  = 


W  B. 


W  B. 


a  = 


AW  B 


AW  B 


a  =  Z  -  A,V 


The  matrices  Ai,  W,  Z,  and  V  are  defined  in  Eqs.  (4.2.6).  The  superscript  ( )-*  indicates  the 
inverse  of  the  complex  conjugate  transpose  matrix.  Equation  (4.2. 10)  can  be  solved  for  X 
by  using  the  techniques  described  in  the  papers  of  Potter  (1966)  or  Martensson  (1971).  In 
general,  there  exist  multiple  solutions  (X's)  to  this  generalized  algebraic  Riccati  equation. 
With  all  solutions  computed,  the  next  step  is  to  decompose  these  solutions  into  Go's  and  Gi  's. 
It  is  shown  in  the  paper  by  Zimmerman  and  Widengren  (1989)  that  for  a  given  solution  X, 
any  selection  of  Gq  and  Gi  satisfying  X  =  Gf  'Gq  results  in  the  same  adjusted  damping 
(Da)  and  stiffness  (Ka)  matrices.  Hence,  either  Gi  (or  Gq)  can  be  chosen  arbitrarily,  as  long 
as  its  inverse  exists.  Then,  Gq  (or  Gi)  is  calculated  from  the  relationship  X  =  Gf  'Gq. 

For  each  set  (Gq,  Gi),  a  feedback  gain  matrix  F  is  calculated  from  Eq.  (4.2.6),  and  the 
corresponding  adjusted  damping  (Da)  and  stiffness  (Ka)  matrices  are  computed  using  Eqs. 
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(4.2.5).  At  this  point,  a  rationale  is  proposed  to  choose  the  most  meaningful  adjusted 
damping  and  stiffness  matrices.  Among  all  computed  sets  (Da,  Kg),  it  is  apparent  that  only 
the  ones  that  are  real  and  symmetric  are  acceptable.  When  dealing  with  a  model  refinement 
problem,  among  all  acceptable  solutions,  the  final  selection  could  be  made  by  choosing  the 
set  (Da,  Ka)  that  minimizes  the  cost  function, 

J  =  q  II  K  -  Ka  If  +    II  D  -  Da  IIf  (4.2.11) 

II  D  If 
where      q  =  t; — ^ 

II  '^IIf 

Clearly,  this  process  selects  the  set  (Da,  Kg)  that  results  in  a  minimum  change  from  the 
original  set  (D,  K).  The  scale  factor  q  in  Eq.  (4.2.11)  is  used  to  give  equal  weight  to  the 
changes  in  D  and  K.  For  the  damage  detection  problem,  there  is  no  unique  rationale  to 
choose  the  "best"  set  (Da,  Ka).  A  physically  intuitive  approach  is  to  use  engineering 
judgement  in  selecting  the  "best"  updated  model.  Thus,  all  acceptable  "adjusted"  sets  of 
solutions  should  be  inspected  to  determine  which  best  provides  information  concerning  the 
state  of  damage. 

4.2.3     Best  Achievable  Eigenvectors 

From  standard  eigenstructure  assignment  theory  (Andry,  et  al.  1983),  it  is  shown  that  the 
measured  eigenvectors  are  not  always  exactly  assignable  to  the  adjusted  finite  element 
model.  In  fact,  it  can  be  shown  that  the  measured  eigenvectors  are  assigned  exacdy  if  and 
only  if  they  lie  in  their  respective  achievable  subspace.  The  achievable  subspace  associated 
to  the  i'*^  mode  is  defined  by 

Lj  =  (Mkf  +  DXj  +  K)-iBo  (4.2.12) 

where  X^  is  the  measured  eigenvalue  of  the  i^^  mode.    When  all  n  components  of  the 
experimental  eigenvectors  are  available,  the  i^^  best  achievable  eigenvectors  is  defined  as  the 
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least  square  projection  of  the  i^^  experimental  eigenvector  Vg.  on  the  i^''  achievable  subspace 
Lj.  This  projection  is  schematically  illustrated  in  Figure  4.1. 


Achievable  Subspace 


Figure  4. 1    Best  Achievable  Eigenvector  Projection. 


This  best  achievable  eigenvector  is  given  by 


Xea,  =  Li 


L,      L: 


LjVp 


(4.2.13) 


When  only  a  subset  s  of  the  eigenvector  components  are  measured,  s  <  n,  the  least  square 
projection  discussed  above  can  be  used  to  simultaneously  expand  and  project  the  measured 
eigenvectors.  In  this  case,  the  i*  expanded  best  achievable  experimental  eigenvector  is 
given  by 


— ea:         ^i 


Li  Lj 


-1 


LiXe, 


(4.2.14) 


where  Lj  are  the  rows  of  Lj  which  correspond  to  the  measured  eigenvector  components. 
Notice  that  the  calculation  of  the  p  achievable  subspace  using  Eq.  (4.2.13)  requires  p 
inversions  of  an  nxn  matrix.  Although  the  matrix  to  be  inverted  is  typically  banded,  this 
may  present  a  practical  computational  burden  when  dealing  with  large  FEMs.  The  next 
section  discusses  an  approach  that  does  not  require  the  actual  computation  of  the  achievable 
subspaces  and  hence  avoids  this  computational  burden. 


■^ai 
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4.2.4     Selection  of  Bn  :  The  Subspace  Rotation  Method 


So  far,  the  control  influence  matrix  has  not  yet  been  completely  defined.  The  preceding 
formulation  suggest  that  different  Bq  may  possibly  result  in  different  adjusted  FEM.  Hence, 
it  is  essential  to  develop  a  physically  meaningful  rationale  to  select  Bq. 

Zimmerman  and  Widengren  (1989,  1990)  proposed  an  approach,  termed  the  mode 
selection  method,  that  consists  in  selecting  Bq  such  that  the  unmeasured  modes  of  the 
structure  are  nearly  unchanged.  In  other  word,  Bq  is  selected  such  that  only  the  measured 
modes  of  the  structure  are  corrected.  This  selection  technique  fixes  the  achievable  subspaces 
in  which  the  eigenvectors  must  lie,  and  hence  places  a  limitation  on  the  assignment  process. 
In  most  studied  cases,  the  experimental  eigenvectors  were  not  assigned  exactly  since  their 
assignment  "success"  depends  on  the  locations  vis-a-vis  the  achievable  subspaces  set  by  the 
selection  of  Bq. 

In  this  work,  a  new  method  of  selecting  Bq,  termed  the  subspace  rotation  method,  is 
proposed.  The  subspace  rotation  method  is  based  on  selecting  Bq  such  that  the  measured 
eigenvectors  lie  exactly  in  the  achievable  eigenvectors  subspaces.  This  procedure  is 
illustrated  in  Figure  4.2  and  is  accomplished  by  setting  Bq  as 

(4.2.15) 


j  =  l, 


where  Vg.  is  the  eigenvector  associated  with  the  j*  experimentally  measured  mode,  and  it 
is  assumed  that  all  n  components  of  the  experimental  eigenvectors  are  available.  This  could 
be  accomplished  by  any  of  the  procedures  discussed  in  Chapter  2.  Clearly,  when  Bq  is 
selected  as  shown  in  Eq.  (4.2.15),  the  measured  expanded  eigenvectors  lie  exactly  in  the 
achievable  subspaces  defined  in  Eq.  (4.2.12).  Hence,  there  is  no  need  for  the  projection 
operations  defined  in  Eq.  (4.2.13).  This  eliminates  the  required  p  inverses  of  nxn  matrices 
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involved  in  computing  the  achievable  subspaces.  As  will  be  seen  in  Chapters  5  and  6,  the 
elements  of  Bq,  as  defined  by  Eq.  (4.2. 15),  give  an  indication  to  the  pseudo-controller  about 
the  extent  of  modification  of  each  DOF  in  order  for  the  structure  to  exhibit  the  j^^  measured 
eigenvalue  and  eigenvector. 


Rotated  Subspace 

■  Nominal  Subspace 


Figure  4.2    Rotation  of  the  Achievable  Subspace. 


4-3    Numerical  Illustrations 

In  this  section,  the  characteristics  of  the  proposed  enhancement  to  the  symmetric 
eigenstructure  assignment  model  refinement  algorithm  (SEAMRA)  are  evaluated  and 
compared  to  other  refinement  techniques  for  two  example  problems.  The  first  problem  is  a 
widely-used  spring-mass  computer  simulated  example  (Kabe,  1985).  It  is  used  here  for  the 
purposes  of  illustrating  model  refinement  for  a  large  local  discrepancy,  analogous  to  a 
damage  detecdon  situation.  The  phenomena  of  global/local  mode  switching  and  load  path 
preservation  are  examined  in  this  problem.  The  second  problem  is  used  to  illustrate  the 
characteristics  of  the  enhanced  SEAMRA  in  updating  the  finite  element  model  (FEM)  of  a 
laboratory  cantilever  beam  using  actual  measured  modal  parameters. 
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4.3.1     Damage  Detection:  Kabe's  Prohlem 

Kabe's  eight  degree  of  freedom  spring-mass  system  is  shown  in  Figure  4.3.  The  mass 
and  stiffness  properties  of  the  system  are  included  in  the  figure.  This  problem  presents  a 
challenging  situation  for  damage  detection  in  that  stiffness  values  of  various  magnitudes  are 
included.  The  model  exhibits  closely-spaced  frequencies  and  both  local  and  global  modes  of 
vibration. 
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Figure  4.3    Kabe's  Problem. 


A  variation  of  Kabe's  original  problem  is  used  here.  Rather  than  the  standard  initial 
model  commonly  used,  which  has  incorrect  values  for  all  of  the  connecting  springs,  only  a 
single  spring  constant  is  changed.  This  is  reflective  of  the  fact  that  damage  may  occur  as  a 
large  local  change  in  the  stiffness  of  a  structural  member. 
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4.3.1.1    Local  to  Global  Mode  Change 

V" 

In  the  first  problem,  Kabe's  initial  model  is  only  incorrect  for  the  spring  between 
masses  3  and  5.  A  value  of  500,  five  times  that  of  the  exact  spring,  is  assumed  in  this  problem. 
Changing  the  spring  value  from  500  to  100  also  causes  a  local  mode  of  vibration  to  be 
replaced  by  a  global  mode,  thus  presenting  a  difficult  challenge  for  damage  detection. 

Figure  4.4  presents  element-by-element  stiffness  matrix  results  for  applying  the  Baruch 
and  Bar  Itzhack  update  (1978)  and  the  symmetric  eigenstructure  assignment  model 
refinement  algorithm.  Baruch  Damage  indicates  that  the  update  was  made  using  Baruch  and 
Bar  Itzhack's  algorithm.  SEA-M  Damage  indicates  that  the  update  was  made  using  the 
SEAMRA  with  Bq  selected  by  using  the  modal  (M)  selection  method.  SEA-SR  Damage 
indicates  that  the  update  was  made  by  using  the  SEAMRA  with  Bq  selected  using  Subspace 
Rotation  (SR)  method.  The  x-coordinate  on  all  plots  are  the  indices  of  a  column  vector 
constructed  by  storing  the  upper  triangular  portion  of  the  stiffness  matrix  in  a  column  vector. 
The  y-coordinate  on  each  plots  consists  of  the  difference  between  the  updated  stiffness 
matrix  elements  and  the  original  stiffness  matrix. 

In  the  first  case,  as  shown  in  Figure  4.4,  it  is  assumed  that  only  the  fundamental  mode  of 
vibration  is  measured,  but  all  eigenvector  components  have  been  measured.  Thus,  no 
expansion  of  eigenvectors  is  required.  It  is  evident  from  Figure  4.4  that  the  Baruch  update  is 
unable  to  discern  the  damage,  but  that  both  the  SEA-M  and  SEA-SR  are  able  to  clearly  locate 
the  damage.  In  fact,  the  SEA-SR  was  able  to  exactly  reproduce  the  correct  stiffness  matrix. 
This  was  true  independent  of  which  mode  was  used  in  the  update.  Also,  it  should  be  noted 
that  the  Baruch  update  tends  to  focus  elemental  changes  in  the  third  and  fifth  row  of  the 
stiffness  matrix,  indicating  the  possibihty  of  damage  between  these  degrees  of  freedom,  but 
certainly  giving  no  clear  indication  to  the  extent  of  damage.  As  is  evident  from  the  plot,  the 
Baruch  update  has  spread  errors  over  several  elements.  Using  the  algorithm  of  Lin  (1990), 
the  damage  vector  is  given  as  a  =  [1.0  0.93  0.72  0.83  0.70  0.90  0.97  1.0]^  where  the 
element  number  corresponds  to  the  structural  DOFs  and  a  number  less  than  1  indicates  the 
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Figure  4.4   Results  for  Kabe's  Problem  using  the  1st  Mode,  Full  Eigenvector. 


possibility  of  damage  affecting  that  DOF.  It  is  obvious  that  DOFs  3  and  5  are  affected  by 
damage,  but  the  results  also  indicate  strong  damage  of  DOF  4. 

In  the  second  case,  as  shown  in  Figure  4.5,  it  is  assumed  that  the  first  three  modes  of 
vibration  have  been  measured,  but  only  the  first  three  components  of  the  eigenvectors  have 
been  measured.  The  eigenvectors  components  were  expanded  for  the  Baruch  update  using 
dynamic  expansion  (Berman  and  Nagy,  1983)  with  subsequent  orthogonalization  (Baruch 
and  Bar  Itzhack,  1978).  The  least  squares  expansion  was  used  for  the  SEA-M  update.  The 
SEA-SR  update  utilized  the  orthogonal  Procrustes  expansion  (Chapter  2).  In  comparing 
Figure  4.5  to  Figure  4.4,  it  is  clear  that  the  damage  detection  capabilities  of  all  three 
algorithms  have  been  degraded  when  using  expanded  mode  shapes,  even  though  more 
modes  have  been  measured.  However,  both  the  SEA-M  and  SEA-SR  updates  give  a  clear 
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Figure  4.5    Results  for  Kabe's  Problem  using  Modes  1,  2,  3, 
and  Eigenvectors  Components  1,  2,  3. 


indication  to  both  the  location  and  extent  of  damage.  Using  Lin's  algorithm,  the  damage 
vector  is  given  as  a  =  [1.0  0.81  0.75  0.83  0.82  0.79  0.85  l.OjT.  It  is  difficult  from 
inspection  of  a  to  determine  the  location  of  damage. 

The  effect  of  applying  the  iterative  load  path  algorithm  described  in  Chapter  2  in  the 
update  procedure  is  shown  in  Figure  4.6.  For  the  Baruch  update,  100  iterations  were 
performed.  For  the  SEA-M  and  SEA-SR  updates  2  and  3  iterations  respectively,  were 
performed.  The  iterations  were  halted  early  for  both  SEA  updates  because  the  discrepancy 
between  the  eigenstructure  before  and  after  masking  was  within  the  numerical  precision  of 
the  symmetric  eigenstructure  assignment  software.  It  is  seen  that  the  load  path  enforcement 
further  enhances  the  damage  detection  capability  of  both  SEA  updates. 
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Figure  4.6    Results  for  Kabe's  Problem  using  Load  Path  Preservation, 
Modes  1,  2,  3,  and  Eigenvectors  Components  1,  2,  3. 


4.3.1.2    Consistent  Modes 


In  the  second  problem,  the  initial  model  is  only  incorrect  for  the  spring  between  masses  4 
and  6.  A  value  of  200,  two  times  that  of  the  exact  spring,  is  assumed  in  this  problem.  In  this 
problem,  all  global  and  local  modes  remain  global  and  local  modes  respectively  after 
damage.  It  should  be  noted  that  finding  a  problem  with  this  feature  was  difficult. 

In  the  first  case,  as  shown  in  Figure  4.7,  it  is  assumed  that  only  the  fundamental  mode  of 
vibration  is  measured,  but  all  eigenvector  components  have  been  measured.  It  is  evident 
from  Figure  4.7  that  the  Baruch  and  SEA-M  update  are  unable  to  discern  the  damage,  but 
that  the  SEA-SR  is  able  to  clearly  locate  damage.  In  fact,  the  SEA-SR  was  able  to  exactly 
reproduce  the  correct  stiffness  matrix.  Again,  this  was  true  independent  of  which  mode  was 
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Figure  4.7    Results  for  Kabe's  Problem  using  the  1st  Mode,  Full  Eigenvectors. 


used  in  the  update.  It  should  be  noted  that  the  Baruch  and  SEA-M  update  tends  to  focus 
elemental  changes  in  the  fourth  and  sixth  rows  of  the  stiffness  matrix,  indicating  the 
possibility  of  damage  between  these  degrees  of  freedom,  but  certainly  giving  no  clear 
indication  to  the  extent  of  damage.  Using  Lin's  algorithm,  the  damage  vector  is  given  as  a  = 
[1.0  0.98  0.92  0.85  0.87  0.84  0.95  1.0]^.  This  algorithm  does  not  clearly  identify  the 
damage  location. 

In  the  second  case,  as  shown  in  Figure  4.8,  it  is  assumed  that  the  first  three  modes  of 
vibration  have  been  measured,  but  only  DOFs  1,3,  and  6  of  the  eigenvectors  have  been 
measured.  In  comparing  Figure  4.8  to  Figure  4.7,  it  is  clear  that  the  damage  detection 
capability  of  all  three  algorithms  has  again  been  degraded  when  using  expanded  mode 
shapes.  Only  the  SEA-SR  update  gives  a  clear  indication  to  the  location  of  damage,  but  is 
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unable  to  predict  the  exact  extent.  Using  Lin's  algorithm,  the  damage  vector  is  given  as  a.=     — 
[0.99  0.18  0.55  0.34  0.52  0.33  0.41   LOf.  Again,  it  is  difficult  from  inspection  of  a  to    ;• 
determine  the  location  of  damage.  In  fact,  inspection  of  a  indicates  that  DOF  2  is  the  most 
likely  damaged  DOF.  «       -■ 
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Figure  4.8    Results  for  Kabe's  Problem  using  Modes  1,  2,  3, 
and  Eigenvectors  Components  1,  2,  6. 


It  should  be  noted  that  in  this  problem,  it  was  critical  to  have  the  proper  DOFs  measured. 
When  the  second  test  case  was  run  with  the  first  three  DOFs  measured,  no  algorithm  was  able 
to  locate  damage.  In  this  case,  the  eigenvectors  components  were  relatively  unaffected  by 
damage,  thus  causing  substantial  error  in  the  eigenvector  expansion  process. 
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4.3.2     Model  Refinement  of  a  Cantilever  Beam:  Experimental  Study 

4.3.2.1    Modal  Test  Description 

The  structure  used  in  this  investigation  is  the  aluminum  cantilevered  beam  shown  in 
Figure  4.9.  The  dimensions  and  material  properties  of  the  beam  are  given  in  Table  4.1. 
Experimental  modal  analysis  of  the  beam  was  performed  to  measure  its  modal  properties. 
Six  equally  spaced  translational  degrees  of  freedom  shown  in  Figure  4.9  were  selected  as 
measurement  locations.  The  modal  properties  of  the  first  three  modes  of  vibration  were 
determined  using  frequency  domain  techniques  and  single  degree  of  freedom  curve  fitting 
algorithms.  The  excitation  source  used  in  the  testing  was  an  impact  hammer  and  the  driving 
point  measurement  was  an  accelerometer  mounted  at  the  free  end  of  the  beam.  Impact  and 
exponential  windows  were  utilized  to  improve  frequency  response  calculations.  At  each 
measured  degree  of  freedom,  five  frequency  response  measurements  were  made  and 
averaged  to  reduce  the  effects  of  measurement  noise.  Natural  frequencies,  damping  ratios 
and  mode  shapes  of  the  beam's  first  three  modes  of  vibration  were  identified  and  are  reported 
in  Tables  4.2  and  4.3. 


1 / 


Figure  4.9    Experimental  Cantilever  Beam. 
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Table  4. 1    Structural  Properties  of  the  Cantilever  Beam. 


Length  -  0.84  m 


Mass/Length  -  2.364  kg/m 


Moment  of  Inertia  -  3.02e-9  m"^ 


Youngs  Modulus  -  70  GPa 


Table  4.2    Measured  Natural  Frequencies  and  Damping  Ratios  of  the 
Cantilever  Beam. 


Mode# 

Natural  Frequency  (Hz) 

Damping  Ratio 

(%) 

1 

7.25 

4.41 

2 

45.55 

0.68 

3 

127.01 

0.33 

Table  4.3    Measured  Mode  Shapes  of  the  Cantilever  Beam. 


Mode 

1 

2 

3 

DOF 

1 

1.00 

1.00 

1.00 

2 

0.95 

0.16 

-0.30 

3 

0.65 

-0.53 

-0.61 

4 

0.36 

-0.72 

0.20 

5 

0.15 

-0.52 

0.75 

6 

0.03 

-0.12 

0.28 

4.3.2.2    Finite  Element  Model  Description 

A  twelve  DOF  undamped  finite  element  model  (FEM)  of  the  beam  was  generated  using 
six  equal  length  beam  elements  as  shown  in  Figure  4.9.  The  beam  element  has  two  degrees  of 
freedom  (DOFs)  at  each  node;  bending  and  rotation.  This  model  was  then  reduced  using 


57 


Guyan  reduction  (1965)  eliminating  the  rotational  degrees  of  freedom.  There  are  several 
possible  errors  affecting  the  accuracy  of  this  FEM.  The  most  obvious  is  the  fact  that  a  perfect 
cantilever  condition  is  assumed.  In  addition,  an  artificial  error  was  purposely  introduced  by 
selecting  the  Young's  Modulus  higher  than  that  often  assumed  for  aluminum. 

4.3.2.3    Application  of  the  SEAMRA. 

Because  the  "true"  finite  element  model  is  unknown,  a  comparison  between  the  "true" 
and  updated  structural  matrices  is  not  possible.  Besides  such  comparison,  a  fair  and  useful 
judgement  on  the  quality  of  SEAMRA  updating  capability  can  be  obtained  by  comparing 
actual  experimental  frequency  response  functions  with  those  predicted  by  the  initial  and 
updated  FEMs.  Figure  4.10  shows  a  comparison  of  frequency  response  functions  measured 
between  degrees  of  freedom  1  and  3  (i.e.  sensor  measurement  at  DOF 1  and  impact  excitation 
at  DOF  3).  The  solid  curve  corresponds  to  the  experimental  data,  while  the  dotted  line 
corresponds  to  that  predicted  by  the  original  analytical  FEM.  It  is  apparent  that  the 
discrepancies  between  the  frequency  response  function  increases  as  the  frequency  of 
excitation  increases.  This  is  in  part  due  to  the  fact  that  the  assumption  of  a  perfect  cantilever 
condition  affects  the  higher  modes  of  vibration  to  a  greater  extent.  The  dashed  lines  in  this 
figure  corresponds  to  the  SEA-SR  updated  finite  element  model.  It  is  clear  from  this 
comparison  that  the  SEA-SR  provided  a  great  deal  of  improvement  to  the  original  analytical 
FEM.  Inspection  of  the  updated  stiffness  matrix  indicates  that  changes  occur  throughout  the 
matrix,  indicating  that  the  discrepancy  between  the  original  and  refined  FEMs  was  due  to 
degradation  of  some  global  structural  property  (Youngs  Modulus),  as  opposed  to  some  form 
of  local  damage,  as  seen  in  the  previous  problem. 

4.4    Discussion  of  the  SEAMRAs  Formulation 

In  some  problems  SEAMRA  in  conjunction  with  either  the  subspace  rotation  or  the 
modal  selection  method,  failed  to  find  a  symmetric  updated  FEM  (symmetric  Da  and  Ka). 
This  shortcoming  was  especially  encountered  in  practical  situations  when  the  experimental 
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Figure  4. 10    Experimental  and  Analytical  Frequency  Response 
Functions  of  the  Cantilever  Beam. 


modal  properties  were  corrupted  by  measurement  errors.  This  can  be  attributed  to  the  fact 
that  there  is  no  symmetric  updated  FEM  that  is  consistent  with  the  present  SEAMRA's 
formulation.  Recall  that  SEAMRA's  modifications  to  the  initial  stiffness  and  damping 
matrices  in  its  updating  process  are  given  by 

AK  =  BoHqBJ         with      Ho  =  FGq 
AD  =  BoHiBJ         with      Hj  =  FGj 


(4.5.1) 


Clearly,  from  Eq.  (4.5.1),  the  perturbations  (i.e.  modifications)  to  the  initial  stiffness  and 
damping  matrices  are  constrained  by  the  relationship 

range(AK)    =  range(AK'^)    =  range(AD)   =  range(AD''')  (4.5.2) 
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This  relationship  can  always  be  satisfied  for  the  cases  when  the  experimental  modal 
properties  are  consistent  with  an  update  in  which  (i)  either  AK  or  AD  is  zero  or  (ii)  AK  is 
proportional  to  AD  (AK  =  p  AD,  p  is  a  scalar).  For  all  other  cases,  SEAMRA  might  fail  to 
produce  symmetric  AK  and  AD.  A  more  flexible  and  general  formulation  that  accounts  for 
such  shortcoming  is  proposed  in  Section  6.5  of  Chapter  6.  The  formulation  as  presented  in 
Chapter  6  is  more  elegant,  efficient  and  guarantees,  for  all  situations,  a  symmetric  updated 
FEM. 

4.5    Summary 

A  previously  developed  model  refinement  algorithm  based  on  the  general  mathematical 
framework  of  eigenstructure  assignment  theory  (Zimmerman  and  Widengren,  1989,  1990) 
has  been  extended  and  improved.  A  technique  to  enhance  eigenvector  assignability  of  the 
algorithm  has  been  developed.  The  method  consists  of  rotating  the  achievable  eigenvector 
subspaces  into  the  experimentally  measured  eigenvectors.  The  subspace  rotation  method, 
used  in  conjunction  with  one  of  eigenvector  expansion  techniques  discussed  in  Chapter  2, 
results  in  both  a  decrease  in  the  computational  burden  as  well  as  an  increase  in  the  accuracy  of 
the  assigned  eigenvectors.  Finally,  the  improved  algorithm  (SEA-SR)  was  tested  for  its 
suitability  for  model  refinement  and  structural  damage  assessment.  The  performances  of 
SEA-SR  in  a  damage  assessment  problem  on  a  challenging  simulated  structure  was 
presented  and  compared  to  other  algorithms.  The  results  acquired  using  the  SEA-SR  were 
superior. 


CHAPTER  5 
DAMAGE  LOCATION:  THE  SUBSPACE  ROTATION  ALGORITHM 


5.1    Introduction 

In  this  chapter,  a  computationally  attractive  algorithm  is  proposed  to  provide  an  insight 
to  the  location  of  structural  damage.  The  proposed  algorithm  is  similar  to  the  Modal  Force 
Error  Criteria  proposed  by  several  researchers  (  Ojalvo  and  Pilon  1988;  Gysin,  1990). 
However,  a  greater  insight  of  the  Modal  Force  Error  criteria  is  provided.  Further,  a  new 
viewpoint  which  allows  for  the  reduction  of  the  effects  of  measurement  errors  in  the 
experimental  modal  parameters  for  a  certain  class  of  structures  is  also  discussed.  As  will  be 
shown  in  the  next  sections,  the  proposed  damage  location  algorithm  requires  only 
matrix-scalar  and  matrix- vector  multiplication. 

5.2    The  Subspace  Rotation  Algorithm:  The  Direct  Method 

Assume  that  an  n-DOF  finite  element  model  of  the  "healthy"  (undamaged)  structure 
exists.  As  seen  in  the  earlier  chapters,  the  standard  differential  equation  governing  the 
dynamic  motion  of  such  structures  is  given  by 

Mw  +  Dw  +  Kw  =  0  (5.2.1) 

where  M,  D,  and  K  are  the  n  x  n  analytical  mass,  damping,  and  stiffness  matrices,  w  is  a  n 
x  1  vector  of  positions  and  the  overdots  represent  differentiation  with  respect  to  time.  The 
eigenvalue  problem  associated  with  Eq.  (5.2.1)  in  second  order  (lambda)  form  is  given  as 

(llM  +  X^p  +  KJv^^  =  0  (5.2.2) 
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where  X^  and  v,^  denote  the  i^^  eigenvalue  and  eigenvector,  respectively,  of  the  pre-damaged 
"healthy"  structure.  It  is  assumed  that  Eq.  (5.2.2)  is  satisfied  for  all  measured  "healthy" 
eigenvalues/eigenvectors.  This  can  be  enforced  by  correlating  the  original  FEM  (M,D,  and 
K),  possibly  through  the  use  of  a  model  refinement  procedure. 

Next,  consider  that  the  p  eigenvalues  and  eigenvectors,  Xji  and  vji,  of  a  post-damage 
modal  survey  of  the  structure  are  available,  in  which  "k^  ^  Xj^ ,  v^  ^  Vj^ .  In  the  present 
formulation,  it  is  assumed  that  the  dimension  of  the  measured  eigenvector  is  the  same  as  the 
analytical  eigenvector.  As  discussed  in  Chapter  2,  this  is  true  (i)  when  all  FEM  DOFs  are 
measured  (ii)  after  the  application  of  an  eigenvector  expansion  algorithm,  or  (iii)  after  the 
application  of  a  finite  element  model  reduction  algorithm.  The  ideal  situation  would  be  to 
measure  all  FEM  DOFs  since  the  eigenvector  expansion  process  would  introduce  additional 
errors  in  the  "expanded"  eigenvectors  and  the  model  reduction  process  would  introduce 
errors  in  the  FEM.  It  should  be  noted  that  in  both  cases  the  additional  errors  may  become 
significant  as  the  ratio  of  measured  to  unmeasured  DOFs  become  smaller. 

Let  AMd,  ADd,  and  AKd  be  the  exact  perturbation  matrices  that  reflect  the  nature  of  the 
structural  damage.  For  localized  damage,  the  exact  perturbation  matrices  are  sparse  matrices 
with  the  nonzero  elements  reflecting  the  state  of  damage.  The  lambda  equation  for  the 
damaged  structure  is,  defined  by 

(X2  (M  -  AMd)  +  X^fp  -  AD^)  +  (K  -  AKd))vd^  =  0         (5.2.3) 
Although  only  p  of  the  n  eigenvalues/eigenvectors  are  assumed  measured,  p  «  n,  Eq.  (5.2.3) 
holds  for  any  particular  eigenvalue  and  eigenvector  of  the  damaged  structure  because  the 
perturbation  matrices  are  assumed  to  be  exact.  Grouping  all  perturbation  matrices  on  the 
right-hand  side  defines  a  damage  vector  dj,  /,. 

d-   =  ZjVj  (5.2.4a) 


=  (X2  AMj  +  Xd_ADd  +  AKd)v^^  (5.2.4b) 


where     Z^    =  XJM  +  X^D  +  K 
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Although  Eqs.  (5.2.4a,  b)  yield  the  same  damage  vector,  it  should  be  noted  that  Z^  and 
the  coefficient  matrix  of  Eq.  (5.2.4b)  are  not  equal.  When  the  measured  eigendata  are  not 
corrupted  by  noise,  an  inspection  of  dj  in  terms  of  the  Eq.  (5.2.4b)  reveals  that  the  j*^  element 
of  di  will  be  zero  when  the  j^^  j-ows  of  the  perturbation  matrices  are  zero,  i.e.  the  finite  element 
model  for  the  j'^^  degree  of  freedom  is  not  directly  affected  by  damage.  Conversely,  a  degree 
of  freedom  whose  finite  element  model  has  been  affected  by  damage  will  result  in  a  nonzero 
entry  in  dj.    Thus,  the  degrees  of  freedom  which  have  been  affected  by  damage  can  be 
determined  by  inspecting  the  elements  of  dj.  Vector  di  as  defined  in  Eq.  (5.2.4b)  also  reveals 
that  only  a  single  mode  of  vibration  needs  to  be  measured  exactly  to  determine  the  damage 
locations  exactly.     This  is  true  in  even  multiple  member  damage  situations.     More 
importantly,  the  vector  dj  can  be  determined  from  the  original  finite  element  model  (M,D,K) 
and  the  measured  eigenvalues  and  eigenvectors,  X^  and  v^,  using  Eq.  (5.2.4a).  Thus,  there  is 
no  need  to  use  a  model  refinement  algorithm  to  attempt  to  estimate  the  exact  perturbation 
matrices  in  order  to  locate  the  damage.    If  the  damping  term  is  ignored,  Eq.  (5.2.4a)  is 
essentially  the  Modal  Force  Error  criteria  as  proposed  by  Ojalvo  (1988)  for  use  as  a 
diagnostic  "tool"  to  locate  modelhng  errors  in  FEMs.    A  physical  interpretation  of  Eq. 
(5.2.4a)  provided  by  Ojalvo  was  "...  di  is  the  applied  harmonic  force  error  distribution, 
applied  at  frequency  X^,  which  is  necessary  to  cause  the  analytical  model  to  vibrate  with 
mode  shape  v^j  ...".  However,  Eq.  (5.2.4b)  provides  a  much  clearer  interpretation  of  the 

damage  vector  di  for  the  damage  location  problem  in  which  the  perturbation  matrices  are 
sparse. 

In  practice,  the  perfect  zero/nonzero  pattern  of  the  damage  vector  di  rarely  occurs  due  to 
errors  present  in  the  experimentally  measured  eigenvalues  and  eigenvectors.  This  scenario 
was  studied  and  discussed  for  the  undamped  case  by  Gysin  (1990)  in  the  context  of 
eigenvector  expansion  techniques.  Gysin  (1990)  observed  that  in  certain  specific  cases  of 
eigenvector  expansion  errors,  the  damage  vector  defined  by  Eq.  (5.2.1a)  may  lead  to 
incorrect  conclusions  concerning  the  location  of  damage.  The  next  section  proposes  a  new 
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viewpoint  which  allows  for  the  reduction  of  the  effects  of  measurement  errors  for  certain 
classes  of  structures. 

5.3    The  Subspace  Rotation  Algorithm:  The  Angle  Perturbation  Method 
In  order  to  provide  an  alternative  view  of  the  state  of  damage,  Eq.  (5.2.4a)  is  rewritten  as 
d{  =4 Yd.   =      114,11  II  ^djl   cos(0i)  (5.3.1) 

where  d]  is  the  j^*'  component  (or  j^*'  DOF)  of  the  i^*^  damage  vector,  z)  is  the  j*  row  of  the 
matrix  Zj   and  6J  is  the  angle  between  the  vectors  z\  and  v^ . 

In  the  case  when  the  measurements  are  free  of  error,  a  zero  d^  corresponds  to  a  Ql  of 
ninety  degrees,  whereas  a  nonzero  d^  corresponds  to  a  6^  different  from  ninety  degrees. 
Errors  in  the  experimental  measurements  of  modal  parameters  will  cause  slight 
perturbations  in  the  angles  0^  that  destroy  the  zero/nonzero  pattern  of  the  damage  vector. 

One  would  initially  expect  that  the  components  of  dj  corresponding  to  the  damaged  DOFs 
would  be  substantially  larger  than  the  other  elements.  However,  by  inspecting  Eq.  (5.3. 1),  a 

large  d,  component  could  be  due  to  a  z)  row  norm  substantially  larger  than  other  rows  of  Z ^ , 

dj  i 

coupled  with  a  slight  deviation  of  6^  from  ninety  degrees  due  to  measurement  noise.  Hence, 
when  dealing  with  a  structure  whose  FEM  results  in   z\  row  norms  of  different  order  of 

d| 

magnitude,  it  is  more  reasonable  to  use  the  deviation  of  the  angles,  61,  from  ninety  degrees 
for  damage  location, 

a)   =  ej   [^]  -  90°  (5.3.2) 

dJ 


where     Q\  =  cos 
1 


-1 


The  angle  6^  is  determined  from  Eq.  (5.3.2)  and  a^  is  the  j*  component  of  ctj. 


1411  II  ^d, 
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5.4    Practical  Issues 

5.4.1     Cumulative  Damage  Location  Vectors 

The  discussion  in  the  previous  section  suggests  that  for  a  given  mode,  the  damage  is 
locatable  if  the  perturbation  to  the  angle,  6J  due  to  the  measurement  error  is  less  than  the 

angle  perturbation  due  to  the  damage.  Hence,  modes  that  are  highly  affected  by  the  damage 
are  expected  to  provide  better  assessment  to  the  location  of  the  damage  when  errors  are 
present  in  the  measured  eigendata.  As  will  be  seen  in  example  problems  (Chapter  7),  certain 
modes  are  more  susceptible  to  a  given  state  of  damage  than  others.  This  is  mainly  due  to  the 
fact  that  different  elements  of  the  structure  have  different  levels  of  contribution  to  the  total 
strain  energy  of  a  given  mode  (Kashangaki  1992).  Furthermore,  a  particular  part  of  the 
structure  usually  has  different  orders  of  strain  energy  contribution  for  different  modes. 
Usually,  if  the  damage  occurs  in  a  region  of  high  strain  energy  for  a  given  mode,  that  mode 
would  be  highly  susceptible  to  the  damage  and,  hence  would  reflect  the  state  of  damage.  To 
accommodate  this  type  of  problem,  when  the  number  of  measured  modes  p  is  greater  than 
one,  two  different  composite  damage  vectors  may  be  defined  as 


a=pX|ail  (5.4.2) 


:1 

In  Eq.  (5.4.1),  the  damage  vectors,  dj,  are  normalized  with  respect  to  their  corresponding 
eigenvectors,  v ^  .  The  reason  for  this  normalization  is  such  that  the  composite  vector  ignores 

the  inherent  "weighting"  of  ||  v^  ||,  which  is  usually  of  different  orders  of  magnitudes  for 
different  measured  modes.  It  should  be  noted  that  in  the  multi-mode  measurement  case,  Eq. 
(5.4.2)  is  preferable  when  the  values  of  ||  z)  \\  are  of  different  orders  of  magnitude  for 

different  measured  modes.  Again,  in  practice,  the  DOFs  affected  by  the  damage  are  expected 
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to  have  substantially  larger  d  or  a.  Finally,  the  damaged  areas  of  the  structure  can  then  be 
located  using  the  knowledge  of  the  "damaged"  DOFs  and  the  connectivity  of  the  FEM. 

It  is  interesting  to  note  that  Eqs.  (5.2.4a,  b)  reveal  an  interesting  relationship  between 
various  model  refinement  algorithms.  Model  refinement  techniques  attempt  to  approximate 
the  exact  perturbation  matrices  by  using  limited  modal  data,  but  do  so  in  different  manners. 
However,  Eq.  (5.2.4a,  b)  indicates  that  if  the  model  refinement  technique  has  satisfaction  of 
the  eigenproblem  as  an  equality  constraint,  the  calculated  perturbation  matrices  AM,  AD, 
and  AK  are  constrained  to  be  related  to  the  original  finite  element  model  M,  D  and  K  and  the 
measured  eigendata  by  Eqs.  (5.2.4a,  b). 

5.4.2     Eigenvector  Filtering  Algorithm 

In  a  modal  survey,  the  errors  associated  with  the  measured  eigenvectors  are  typically 
greater  than  the  error  associated  with  the  measured  eigenvalues.  In  addition,  in  the  case  of 
incomplete  eigenvector  component  measurements,  these  measurement  errors  are  often 
compounded  with  eigenvector  expansion  induced  errors.  A  simple  eigenvector  noise 
filtering  algorithm  is  proposed  assuming  the  measured  eigenvalues  to  be  correct.  From  the 
cumulative  damage  vector  defined  in  Eqs.  (5.2.4a)  or  (5.2.4b),  and  the  original  FEM 
connectivity,  the  engineer  can  deduce  which  DOFs  have  been  damaged.  It  is  reasonable  to 
assume  that  nonzero  elements  in  each  damage  vector  di  associated  with  "undamaged  DOFs" 
are  due  to  eigenvector  errors.  These  elements  can  then  be  set  to  zero.  In  addition,  the 
magnitude  of  the  elements  of  dj  at  the  "damaged"  DOFs  can  be  adjusted  by  using  knowledge 
of  DOE  connectivity  and  the  properties  of  the  element  property  matrices  connecting  the 
"damaged  DOFs."  The  element  property  matrices  provide  constraints  relating  the  effect  of 
damage  on  each  element  DOF.  The  noise  filtering  algorithm  consists  simply  of  replacing  the 
dj    vectors  by  dr,  where  dr  is  obtained  from  d;  as  described  above.    The  i^^  filtered 

i  ^i 

eigenvector,  v^^,  can  then  be  obtained  from  solving 


(XJM  +  XjD  +  KJvjf^   =  df_  (5.4.3) 
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using  Gaussian  elimination.  In  this  calculation,  the  handedness  of  typical  FEM  matrices 
should  be  exploited.  Essentially,  the  filtered  eigenvector  is  just  the  eigenvector  that  if 
measured  would  have  produced  the  damage  vector  d^ .  Experience  gained  in  using  the 
eigenvector  filtering  algorithm  indicates  that  it  is  best  to  use  structural  matrix  properties  (M, 
D,  K)  that  (i)  are  finite  element  consistent,  and  (ii)  have  not  been  "corrupted"  by 
measurement  noise.  By  finite  element  consistent,  it  is  meant  that  the  property  matrices  can 
be  achieved  by  a  finite  element  program.  Note  that  measurement  noise  can  be  introduced 
in  the  property  matrix  through  a  FEM  refinement  algorithm.  Hence,  the  property  matrices 
that  should  be  used  are  the  original  property  matrices  (unrefined).  The  present  eigenvector 
filtering  algorithm  can  be  useful  in  improving  the  damage  extent  assessment.  This  concept 
<■      is  discussed  in  Chapters  6  and  7.  ,   '        ,     •  \: 

5.5    Summary 

A  computationally  attractive  algorithm  to  determine  the  location  of  damage  in 
structures  was  developed.  The  algorithm  completely  bypasses  the  general  framework  of  the 
model  refinement  problem  and  involves  only  matrix-scalar  and  matrix-vector 
multiplicafions.  The  effect  of  measurement  error  in  the  eigendata  was  discussed  and 
techniques  to  reduce  these  effect  were  presented.  Furthermore,  a  simple  eigenvector 
filtering  algorithm  was  developed.  Practical  example  problems  to  illustrate  and  evaluate  the 
performance  of  the  developed  algorithm  will  be  presented  in  Chapter  7. 
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CHAPTER  6 
THE  MINIMUM  RANK  PERTURBATION  THEORY 


6.1    Background 

The  theory  developed  in  Chapter  5  is  limited  to  determining  the  location  of  structural 
damage.  In  a  practical  situation,  it  is  essential  to  determine  the  extent  of  the  damage  to  get  a 
good  estimate  about  the  overall  integrity  of  the  structure.  In  general,  the  extent  problem,  as 
discussed  in  Chapter  1,  has  been  approached  by  several  researchers  using  existing  model 
refinement  algorithms.  The  formulations  of  these  algorithms  were  obviously  based  on  the 
"model  refinement  philosophy":  minimum  change  made  to  the  original  FEM.  The  minimum 
change  constraint  has  a  clear  tendency  to  smear  the  changes  throughout  the  entire  FEM. 
However,  in  most  cases,  this  philosophy  is  not  consistent  with  the  effect  of  structural  damage 
on  FEMs.   In  fact,  the  effects  of  structural  damage  on  FEMs  are  usually  "non-minimal" 
localized  perturbations.  Structural  damage  often  occurs  at  discrete  locations.  The  effect  of 
damage  on  the  analytical  model  is  often  restricted  to  just  a  few  elements  of  the  finite  element 
model.   The  rank  of  each  element  mass,  damping  or  stiffness  matrix  is  dependent  on  the 
number  of  degrees  of  freedom  defined  by  the  element  and  the  shape  functions  utilized. 
However,  it  should  be  noted  that  in  general  the  element  matrices  are  not  of  full  rank.  For 
example,  the  rank  of  the  6x6  element  stiffness  matrix  of  a  three  dimensional  truss  element  is 
just  one.  Thus,  instead  of  using  the  matrix  Frobenius  norm  minimization  formulation  to 
arrive  at  unique  perturbation  matrices,  minimum  rank  perturbation  constraints  are  enforced. 
In  this  Chapter,  a  computationally  attractive  damage  extent  algorithm  is  proposed.  The 
proposed  damage  extent  algorithm  is  a  minimum  rank  perturbation,  which  is  consistent  with 
the  effects  of  many  classes  of  structural  damage  on  a  FEM. 
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Assume  that  "p  damaged"  eigenvalues  and  eigenvectors  have  been  measured  and  that 
the  original  FEM  has  been  corrected  such  that  its  modal  properties  match  the  measured 
modal  properties  of  the  healthy  model.  The  eigenvalue  problem  of  a  damaged  structure 
shown  in  Eq.  (5.2.3),  for  all  p  measured  modes,  can  be  written  in  matrix  form,  as 

MVjAj  +  DV,A,  +  KV,  =  AMdV,A2  +  AD^V^A,  +  AK^V^  =  B  (6.1.1) 

where      A^  =  dmg[l^^,k^^,  ...,k^^j 


Vd  = 


B     = 


Vd,,Vd,...,Vd 


di,d2,...,dp 


where  all  variables  have  the  same  definitions  as  in  the  previous  chapter.  Note  that  matrix 
B  can  be  determined  from  the  FEM  (M,  D,  K)  and  the  "p"  measured  eigenvalues  and 
eigenvectors.  As  discussed  earlier,  the  damage  extent  problem  consist  of  finding  the 
perturbation  matrices,  AMj,  ADd,  and  AK^,  such  that  Eq.  (6.1.1)  is  satisfied.  As  already 
discussed  in  Chapter  2,  structures  can  be  modelled  using  either  undamped,  proportionally 
damped  or  nonproportionally  damped  finite  element  models.  The  proposed  extent  algorithm 
is  formulated  to  accommodate  all  three  types  of  structural  models.  For  each  type  of  model, 
several  scenarios  of  damage  effects  are  considered.  Practical  issues  that  can  be  used  to 
improve  the  damage  extent  estimate  are  also  presented. 

6.2    The  Minimum  Rank  Perturbation  Theory:  Theoretical  Background 

In  this  section,  the  theoretical  foundation  of  the  Minimum  Rank  Perturbation  Theory 
(MRPT)  is  derived.  This  theory  will  be  extensively  used  throughout  the  remainder  of  this 
chapter. 

PROPOSITION    6.1    Suppose    that    X,  Y  G  R"'^p   are    given    where    p    <   n       and 

rank(X)=rank(Y)=p. 

Define  %  to  be  the  set  of  matrices  A  in  R"™  that  satisfy. 


...■^ 
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AX  =  Y        with        A^^  =  A  (6.2.1) 

Then, 

(La)  If  the  set  3G  is  nonempty,  the  minimum  rank  of  any  matrix.  A,  in  %  is  p. 

Next,  define  %P  to  be  a  subset  of  %  comprised  of  all  A  such  that  rank(A)  =p. 
Then 

(l.b)  If  the  matrix  Y'^X  is  symmetric,  then  one  member  of  %p  is  given  by 


AP  =  YHY'^  (6.2.2) 


with        H  =  (y^^x) 


1 


and 

(l.c)  The  matrix  defined  by  Eq.  (6.2.2)  is  the  unique  member  of  3GP. 

Proof: 

To  prove  Proposition  (6.1. a),  note  that  Eq.  (6.2.1)  is  exactly  satisfied  if  and  only  if 
range(Y)  is  included  in  range(A),  which  is  also  the  range(A'^)  by  symmetry.  This  implies  that 
rank(Y)  =  p  <  rank(A). 

/// 

To  investigate  Proposition  (6. l.b),  assume  that  the  expanded  singular  value 
decomposition  of  one  member,  AP-J,  of  %P  to  be  of  the  form 

jT 


APJ  =  UJ2JUJ  (6.2.3) 

where  UJ  =  [uJ,  u^,  ...  ,uJ,] 

2J  =  diag(oJ,  a^,  ...  .aj,) 
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where  the  superscript  j  indicates  the  j^^  family  member  of  %P,  the  uJ  are  the  left  and  right 
singular  vectors  and  the  ol  are  the  non-zero  singular  values  of  APJ.  In  the  expanded  singular 
value  decomposition,  the  (p+1)  to  n  singular  vectors  are  not  shown  in  the  factorization 
because  of  their  corresponding  zero  singular  values.  Note  that  the  left  and  right  singular 
vectors  are  the  same  because  AP'J  is  restricted  to  be  symmetric.  For  Eq.  (6.2. 1 )  to  be  satisfied, 
the  range  of  Y,  APJ,  APJ""  and  UJ  must  be  equal.  Thus,  any  column  of  Y  can  be  written  as 
a  linear  combination  of  the  uJ's.  The  matrices  Y  and  UJ  are  then  related  by  a  unique  pxp 
invertible  matrix  Q), 

Y  =  UJ  QJ  (6.2.4) 

Substituting  Eq.  (6.2.4)  into  Eq.  (6.2.3)  gives 

APJ  =  y(qJ"'2JQJ"')yT'  =  YHJyT  (6.2.5) 

Thus,  each  family  member  is  uniquely  defined  by  the  factorization  of  Eq.  (6.2.3).  From  Eq. 
(6.2.5),  it  is  evident  that  HJ  is  of  full  rank  because  its  inverse  exists  (HJ~'  =  QJ^2J"'QJ). 

Inspection  of  Eq.  (6.2.5)  reveals  that  the  only  unknown  term  in  the  factorization  is  HJ. 
By  using  the  factorization  of  APJ  as  defined  by  Eq.  (6.2.5),  Eq.  (6.2.5)  can  be  rewritten  as 

Y  =  APJX  =  (yHJyT)x  =  y(hJYTx)  (6.2.6) 

Equation  (6.2.6)  is  satisfied  if  and  only  if  HJ'y'^X  =  Ipxp,  where  Ipxp  is  the  pxp  identity 
matrix.  This  is  true  because  Y  and  X  are  of  full  column  rank.  Thus,  HJ  is  uniquely  calculated 
to  be 

HJ  =  (y^x)"'  (6.2.7) 

/// 
Proposition  (6. 1  .c)  follows  immediately  by  inspecting  the  right  hand  side  of  Eq.  (6.2.7). 
Inspection  reveals  that  HJ  is  the  same  for  all  members  of  %P.  This  fact,  in  conjunction  with 
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Eq.  (6.2.5)  leads  to  the  conclusion  that  AP'J  is  the  unique  member  of  the  set  %P.  This  member 
is  given  by  Eq.  (6.2.2). 

/// 

■ 

At  this  point,  the  MRPT  as  defined  in  Proposition  6. 1  assumes  that  the  matrices  X  and  Y 

are  of  full  rank.  In  practical  uses,  as  will  be  seen  later,  matrix  X  is  usually  of  full  rank.  The 

rank  requirement  on  matrix  Y  can  be  of  some  concern  since  it  is  directly  related  to  the  rank  of 

matrix  A.  The  next  proposition  addresses  the  case  in  which  matrix  Y  is  rank  deficient. 

PROPOSITION  6.2  Suppose  that  X,  Y  G  M"'^'"  are  given  and  rank(X)=m  and  rank(Y)=p, 
where  p  <  m  <  n.  Further,  suppose  that  the  matrix  Y^X  is  symmetric. 
Define  11  to  be  the  set  of  matrices  A  in  R"""  that  satisfies  the  problem, 

APXP  =  YP  with        (AP)"^  =  AP  (6.2.8) 

where  the  superscript  p  indicates  a  rank  p  matrix.    In  Eq.  (6.2.8),  XP,  YP  G  R"^p  are 

corresponding  full  rank  submatrices  of  X  and  Y. 

Then 

The  set  ^contains  a  single  member.  A?,  that  can  be  calculated  from  Eq.  (6.2.2)  using 
any  corresponding  XP  and  YP. 

Proof: 

The  j*  member  of  the  set  Rl  is  given  by 

AP'J  =  YPJHJYPJ''  (6.2.9) 

with        HJ  =  (ypJ^XPJ)' 

where  the  additional  superscript  (  )-'J  indicates  the  j*  member  of  '^L.  Note  that  HJ  is 
symmetric  since  Y^X  is  symmetric.  The  range  of  any  YP  is  equal  to  the  range  of  Y,  thus 
the  YP'i  and  YP-J  are  related  by 
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YP-'  =  YP'J  Q'J  (6.2.10) 

where  Q*J  G  RP"?  and  rank(Qi'J)  =  p.  By  utilizing  Eq.  (6.2.10),  Eq.  (6.2.9)  can  be  written 
for  the  i^^  member  of  Rl  as 

AP-i     =    YPd(Q'dH*Q*J^)YPJ^  (6.2.11) 

with      H'  =  (q'jV'J^XP'') 

Again,  H'  is  symmetric  because  Y^X  is  symmetric.  In  comparing  Eqs.  (6.2.9)  and  (6.2. 1 1), 
it  is  seen  that  APJ  =  AP''  if 

HJ  =  Q'J  H*  Q'J^  (6.2.12) 

or  equivalently, 

HJ"'  =  Q'J~^H'"^Q'J"'  (6.2.13) 

where  Eq.  (6.2.13)  makes  use  of  the  symmetry  of  H'.  By  using  the  definitions  of  H"  and  HJ, 
Eq.  (6.2.13)  can  be  written  as 

ypJ'^XPJ   =  Q>'J'''XP''VP'J  (6.2.14) 

Pre-multiplying  Eq.  (6.2.14)  by  Q'J^  and  utilizing  the  relation  in  Eq.  (6.2.10)  gives  the 
condition  such  that  AP'J  =  AP-',  namely 


yp-'^xP'J  =  xP'^'^ypj  (6.2.15) 

This  condition  is  clearly  satisfied  since  Y^X  is  symmetric. 

The  actual  uses  of  Proposition  6. 1  and  Proposition  6.2  will  be  clearly  seen  in  the  next 
sections.  The  practical  implication  of  Proposition  6.2  is  discussed  in  detail  in  Section  6.6. 

6.3    Damage  Extent:  Undamped  Structures 

In  some  cases,  the  damping  of  the  system  under  consideration  is  assumed  to  be 
negligible.  For  this  type  of  system,  MRPT  based  algorithms  will  be  developed  assuming  that 
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the  structural  damage  affects  (i)  only  the  mass  properties,  or  (ii)  only  the  stiffness  properties, 
or  (iii)  simultaneously  the  mass  and  stiffness  properties. 

6.3.1     Damage  Extent:  Mass  Properties 

In  this  case,  it  is  assumed  that  the  effect  of  damage  on  the  stiffness  properties  of  the 
structure  is  negligible.  With  this  assumption,  Eq.  (6.1.1)  can  be  rewritten  as 

MVjA^  +  KVd  =  AMdVjAj  ^  B  (6.3.1) 

Note  that  the  eigenvectors  are  real  and  the  eigenvalues  are  purely  imaginary.  Further,  the 
eigenvectors  are  linearly  independent,  which  implies  that  the  matrix  product  V^A^  is  of  full 
column  rank  if  rigid  body  modes  are  not  included.  Assume,  for  the  moment,  that  B  is  of  full 
rank  (rank(B)  =  p).  Then,  Proposition  6.1  can  be  applied  to  determine  the  perturbation 
matrix,  AM^,  as 

AM^  =  B(B'^VdA2)  ~  b"^  (6.3.2) 

by  letting  Y=B  and  X=VjAj .  Note  that  the  required  inversion  is  that  of  a  pxp  matrix,  where 
"p"  is  the  number  of  measured  modes.  As  discussed  in  Proposition  6. 1 ,  this  inversion  is 
feasible  if  matrix  B  is  of  full  rank  and  the  rigid  body  modes  of  the  system  are  omitted  in  the 
computations.  When  matrix  B  is  rank  deficient.  Proposition  6.2  should  be  used  to  render 
the  computation  possible. 

The  properties  associated  with  AMj  as  computed  in  Eq.  (6.3.2)  are  as  follows: 

PROPERTY  6.3.1  The  perturbation  matrix,  AMj,  defined  in  Eq.  (6.3.2)  will  be  symmetric  if 
the  eigenvectors,  Vj,  are  stiffness  orthogonal,  i.e.,  the  eigenvectors  are  orthogonal  with 
respect  to  the  original  stiffness  matrix,  K. 
Proof: 

Proposition  (6. 1  .c)  in  conjunction  with  Proposition  (6. 1  .b)  implies  that  the  existence  of 

the  unique  symmetric  rank  p  AMj  requires  the  symmetry  of  the  matrix  product  B^V^jA^. 
The  symmetric  equivalence  associated  with  this  matrix  product  is 
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BTVjA^    =    AjvjB  (6.3.3) 

Substituting  the  expression  for  B  from  Eq.  (6.3.1)  into  Eq.  (6.3.3)  gives 

AjvjMV^Aj  +  VjKV,A2  ^  A^vJMV.A^  +  A^vJKV,  (6.3.4) 

where  the  symmetry  of  M,  K  and  A^  has  been  used  in  writing  Eq.  (6.3.4).  From  Eq.  (6.3.4), 
it  is  clear  that  the  equivalence  is  true  if 


(VJKVJA^    .    A2(VJKVJ  (6.3.5) 


Equation  (6.3.5)  will  obviously  be  satisfied  if  the  measured  eigenvectors  are  stiffness 
orthogonal.  Baruch  (1978)  treated  one  approach  to  mass  orthogonalize  the  measured 
eigenvectors.  A  similar  approach  can  be  used  to  orthogonalize  the  measured  eigenvectors 
with  respect  to  the  stiffness  matrix. 


PROPERTY  6.3.2  The  updated  finite  element  model  (FEM)  defined  by  the  original  mass 
and  stiffness  matrix  along  with  the  perturbation  mass  matrix  computed  using  Eq.  (6.3.2) 
preserve  the  rigid  body  characteristics  of  the  original  FEM. 

Proof: 

This  is  apparent  in  that  the  original  stiffness  matrix  is  unchanged  and  that  the  rigid  body 
modes  are  defined  as  modes  whose  eigenvectors  lie  in  the  null  space  of  the  stiffness  matrix. 

6.3.2    Damage  Extent:  Stiffness  Properties 

Here,  it  is  assumed  that  the  effect  of  damage  on  the  mass  properties  of  the  structure  is 
negligible.  With  this  assumption,  Eq.  (6. 1 . 1)  can  be  rewritten  as 

MV^Al  +  KVj  =  AKjVj  =  B  (6.3.6) 

For  this  problem,  the  eigenvectors  are  real  and  the  eigenvalues  are  purely  imaginary.  The 
eigenvectors  are  also  linearly  independent,  which  implies  that  matrix  Vj  is  of  full  rank,  ff 
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matrix  B  is  assumed  to  be  of  full  rank  (rank(B)=p),  Proposition  6. 1  can  be  used  to  determine 
the  perturbation  to  the  original  stiffness  matrix, 

AK,  =  B(bTv,)"  B^  (6.3.7) 

This  expression  for  AKj  is  determined  by  setting  Y=B  and  X=V(i  in  Eq.  (6.1.3). 
The  properties  associated  with  AKj  as  computed  by  Eq.  (6.3.7)  are  as  follow. 

PROPERTY  6.3.3  The  matrix  AK^  will  be  symmetric  if  the  eigenvectors  are  mass 
orthogonal,  i.e.,  the  eigenvectors  are  orthogonal  with  respect  to  the  original  mass  matrix. 
The  proof  of  Property  6.3.3  follows  very  much  the  same  pattern  as  the  one  presented  for  AMjj 
(Property  6.3.1). 

PROPERTY  6.3.4  The  updated  FEM  defined  by  the  original  mass  and  stiffness  matrices  and 
the  perturbation  stiffness  matrix,  AK^i,  preserves  the  rigid  body  mode  characteristics  if  the 
measured  eigenvectors  and  the  rigid  body  modes  are  mass  orthogonal.  ■ 

Proof: 

The  original  rigid  body  modes  of  an  undamped  system  are  defined  by  the  eigenvalue 
problem, 

KVr  =   XrMVr  =   0  (6.3.8) 

where  the  subscript  r  denotes  the  rigid  body  mode(s)  and  \  is  equal  to  zero.  Thus,  the  rigid 
body  modes  lie  in  the  null  space  of  the  original  stiffness  matrix.  The  rigid  body  modes  of 
the  system  will  be  preserved  in  the  updated  model  if  the  original  rigid  body  modes  lie  in  the 
null  space  of  the  updated  stiffness  matrix, 

e  =  (k  -  AkJv^  (6.3.9a) 

=  -  AKjVr  (6.3.9b) 

where  vector  e  is  zero  if  the  the  rigid  body  modes  are  preserved.  Equation  (6.3.8)  has  been 
used  to  arrive  at  the  expression  shown  in  Eq.  (6.3.9b).  Substituting  Eq.  (6.3.7)  into  (6.3.9b) 
gives 
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e  =  -  B(B'^Vd)     B\,  (6.3.10) 

By  utilizing  the  symmetry  of  the  original  mass  and  stiffness  matrices,  along  with  Eq.  (6.3.6), 
Eq.  (6.3.10)  can  be  expanded  as 

e  =  -  B(bT"V,)"'[vJ"Kv,  +  A.VjMv,]  (6.3.11) 

The  first  term  in  the  parenthesis  is  zero  because  the  matrix-vector  product  Kvr  is  zero  by 
definition.  The  second  term  will  be  zero  if  the  rigid  body  modes  and  the  measured  mode 
shapes  are  mass  orthogonal. 

^^ 

6.3.3     Damage  Extent:  Mass  and  Stiffness  Properties 

In  this  case,  it  is  assumed  that  the  structural  damage  affects  simultaneously  the  mass  and 
stiffness  properties  of  the  structure.  With  this  assumption,  Eq.  (6.1.1)  can  be  rewritten  as 

MV,A2  +  KV,  =  AM.VjAj  +  AKjVd  =  B  (6.3.12) 

6.3.3. 1    Application  of  The  MRPT 

Assume  that  Eq.  (6.1.12)  can  be  decoupled  as  follows, 

AMjVd    =    Bn,  (6.3.13a) 

AK^V^    =    Bk  (6.3.13b) 

Then,  the  MRPT,  as  formulated  in  Proposition  6.1,  can  be  applied  to  determine  the 
perturbation  matrices  AMj  and  AK^j,  as 

AM,  =  B„(bTv,)"'bT  (6.3.14a) 

AKd  =  B,  (bTvO~'bJ  (6.3.14b) 

Note  that  the  matrices  B^Vj  and  BJ^Vj  are  invertible  if  B^  and  B^  are  of  full  rank.  When 
these  rank  requirements  are  not  met.  Proposition  6.2  can  be  used  to  make  the  computations 
possible. 
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6.3.3.2    Decomposition  of  Matrix  B 

The  decomposition  problem  as  illustrated  in  the  previous  section  is  equivalent  to  the 
problem  of  solving  for  the  matrices  Bm  and  B^.  So  far,  the  only  constraint  that  these 
unknown  matrices  must  satisfy  is  given  by  the  expression, 

B    =    Bm  Aj  +  B,,  (6.3.15) 

which  results  from  Eqs.  (6.3. 12)  and  (6.3. 13).  Naturally,  there  is  an  infinite  set  of  solutions 
(Bm,  Bk)  that  satisfy  Eq.  (6.3.15).  To  arrive  at  a  unique  solution,  additional  physically 
meaningful  constraints  can  be  enforced.  The  decomposition  proposed  herein  exploits  the 
cross-orthogonality  relations  that  arise  from  the  symmetric  nature  of  the  property  matrices 
and  the  undamped  assumption.  By  measuring  mass  normalized  "damaged"  eigenvectors 
(which  is  possible  if  a  driving  point  measurement  is  made),  the  cross-orthogonality  relations 
associated  with  the  damaged  structure  can  be  written  as 

vJ(M  -  AMdjVj  =  Ipxp  (6.3.16a) 

VT(K  -  AKJV,  =  diag(a)d,2,    ■  ■  ■  ,0)^^')  =  Q^  (6.3.16b) 

in  which  coj  is  the  natural  frequency  of  the  i^*'  mode  of  the  "damaged"  structure.  Matrix 
Ipxp  is  the  pxp  identity  matrix.  A  rearrangement  of  Eq.  (6.3. 16)  yields 

VjAM^Vd    =    VjMV,    -    Ipxp    =    VjBm  (6-3.17a) 

VjAKjVd    =    VjKV,    -    Q,    =    VjB,  (6.3.17b) 

Clearly,  the  matrices  B^  and  Bk  can  be  computed  from  Eqs.  (6.3.17).  In  the  rare  situation 
that  the  number  of  measured  modes  is  equal  to  the  number  of  DOFs  in  the  FEM  (p  =  n),  these 
can  be  computed  by  simply  inverting  matrix  Vd.  Unfortunately,  as  discussed  earlier,  the 
number  of  measured  modes  is  usually  much  less  than  the  number  of  FEM  DOFs  (p  «  n). 
In  this  case,  the  solution  that  naturally  comes  to  mind  is  to  use  the  pseudo-inverse  of  matrix 
Vj .  The  inconvenience  of  this  approach  is  that  the  sparsity  pattern  of  matrix  B  will  not  be 
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reflected  in  the  computed  matrices  Bm  and  B^.  Remember  that  the  sparsity  pattern  of  B,  as 

discussed  in  Chapter  5,  indicates  the  location  of  the  damage  affecting  the  structure.  A  more 

j-%\.^  physically  intuitive  approach  is  to  constrain  Bm  and  B^  to  exhibit  the  same  sparsity  pattern 

^?  ■  as  matrix  B.  This  is  done  by  casting  B  in  an  equation  similar  to  the  expressions  of  Eqs. 

,  (6.3. 17).  The  problem  in  question  is  then  to  find  an  nxp  matrix  P  that  satisfies 

p(vJb)    =    B  (6.3.18) 

Matrix  P  can  be  computed  as 

P    =    B  (vJb)~  (6.3.19) 


The  inverse  involved  in  this  computation  is  that  of  a  pxp  matrix  which  is  invertible  if  matrix 
B  is  of  full  rank.  Now  that  P  is  computed,  Bm  and  B^  can  be  computed  using  Eq.  (6.3. 19) 
as 


Bn,    =    P(VJMV,    -    Ip,p)  (6.3.20a) 

Bk    =    P(VIKV,    -    qJ  (6.3.20b) 


It  is  clear  from  Eq.  (6.3. 19)  that  P  will  have  the  same  sparsity  pattern  as  matrix  B.  Hence 
Bm  and  Bk  will  also  reflect  the  important  sparsity  pattern  of  B.  The  computed  matrices  Bm 
and  Bk  can  also  be  used  to  determine  the  effect  of  the  damage,  respectively,  on  the  mass,  and 
stiffness  properties.  As  in  Chapter  5,  cumulative  vectors  associated  to  Bm  and  Bk  can  also  be 
defined  when  more  than  one  measured  mode  is  available. 

dk-^y^^  (6.3.21b) 

where  dm  and  d,^  are,  respectively,  the  i"^  column  of  matrix  Bm  and  Bk. 
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PROPERTY  6.3.5  The  perturbation  matrices  (AMd,  AKd)  computed  from  the  MRPT  using 
the  Bm  and  Bk  resulting  from  the  decomposition  discussed  above  will  be  symmetric. 

Proof: 

The  perturbation  matrices  AMj  and  AKj  will  be  symmetric  because  they  satisfy  the 
relationships  in  Eqs.  (6.3.17)  and  the  right  hand  sides  of  these  equations  are  symmetric. 

6.4  Damage  Extent:  Proportionally  Damped  Structures 

Since  many  structures  have  non-negligible  damping,  it  is  of  practical  interest  to  extend 
the  MRPT  to  address  damped  structures.  In  this  analysis,  the  structure  under  consideration  is 
f  assumed  to  exhibit  proportional  damping. 

6.4.1     Damage  Extent:  Stiffness  and  Damping  Properties 

It  is  assumed  that  the  effect  of  the  structural  damage  on  the  mass  properties  is  negligible. 
In  this  context,  Eq.  (6. 1 . 1)  is  rewritten  as 

MV,A2  +  DVjAd  +  KVj  =  AD^VjAd  +  AK^V^  ^  B  (6.4.1) 

The  complex  conjugate  of  Eq.  (6.4. 1)  is 

AD,V,A,  +  AKjVd  =  B  (6.4.2) 

where  the  overbar  indicates  the  complex  conjugate  operator,  and  the  fact  that  ADj,  AK^  and 
Vd  are  real  has  been  used  in  writing  Eq.  (6.4.2).  Subtracting  Eq.  (6.4.2)  from  Eq.  (6.4.1) 
gives 

-        ;  AD,Vd(Ad  -  Aj)  =  (B  -  B)  (6.4.3) 

If  (B  -  B)  is  assumed  to  be  of  full  rank,  Proposition  6. 1  can  be  applied  to  determine  the 
perturbation  matrix,  AD^,  as 

ADd  =  (B  -  B)Hd(B  -  B)"" 


i 


with   Hj  = 


_i  (6.4.4) 


(B  -  B)  V,(A,  -  A,) 
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Note  that  ADj  as  defined  by  Eq.  (6.4.4)  is  real.  Post-multiplying  Eq.  (6.4. 1)  by  A^  and 
Eq.  (6.4.2)  by  A^j  and  subtracting  the  two  equations  leads  to 

AKdV,(Ad  -  A,)  =  (BAd  -  BAj)  (6.4.5) 

where  the  fact  that  A^  and  A^  are  diagonal  matrices  has  been  used  in  writing  Eq.  (6.4.5). 
If  (b Aj  -  B Aj)  is  assumed  to  be  of  full  rank,  Proposition  6. 1  can  also  be  applied  to 
determine  the  perturbation  matrix,  AIQj,  as 


AK,  =  (BA,  -  BA,)H,(BAd  -  BA,f 


with   H|^  = 


(6.4.6) 


"    -       _      T      _  V' 

(BA,  -  BA,)  V,(A,  -  AjJ 

Note  that  AIQ  as  defined  by  Eq.  (33)  is  also  real. 

PROPERTY  6.4. 1  The  perturbation  matrices  ADj  and  AIQ,  as  computed  above,  will  be 
symmetric  if  the  measured  eigenvectors,  Vd,  are  mass  orthogonal;  i.e.,  the  eigenvectors  are 
orthogonal  with  respect  to  the  original  unperturbed  mass  matrix.  i  ~  v^  *'  ' 

Proof: 

Matrix  ADj  is  symmetric  if  Ra  is  symmetric  or,  equivalently,  if  H J '  is  symmetric. 
Hence,  to  get  a  symmetric  ADj,  the  following  equivalence  must  be  satisfied. 

(B  -  BjV.lAd  -  A,)   ^  (A,  -  A,)VJ(B  -  B)  (6.4.7) 

Substituting  the  expressions  for  B  and  B,  from  Eqs.  (6.4.1)  and  (6.4.2) ,  respectively,  into 
Eq.  (6.4.7)  yields 


(aJvJm  +  A,VjD  -  AdVjM  -  A.VJd)  V,(A,  -  A,) 

=  (A,  -  A,)VJ  (mV,A2  +  DV,A,  -  MV,Ad  +  DV,A,) 


(6.4.8) 
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Note  that  in  Eq.  (6.4.8)  the  terms  involving  matrix  AK<i  canceled  out.   A  further  expansion 
and  simplification  of  Eq.  (6.4.8)  yields 

(a2  -  Ad)  VjMV,  (A,  -  A,)   ^  (A,  -  A,)  vJmV,  (a^  -  Ad)    (6.4.9) 


which  is  clearly  satisfied  if  the  measured  "damaged"  eigenvectors,  Vj,  are  mass  orthogonal. 

/// 
Likewise,  the  perturbation  matrix  AK^  as  computed  in  Eq.  (6.4.6)  is  symmetric  if  Hk  is 
symmetric  or,  equivalently,  if  H^ '  is  symmetric.  This  symmetry  requirement  yields  the 
following  equivalence. 

(BA,  -  Ba/v,(A,  -  A  J   =  (A,  -  AjvJ(BA,  -  BAJ        (6.4.10) 
Substitution  of  the  expressions  for  B  and  B  into  Eq.  (6.4.10)  yields 

(AdA^vjM  +  A.VjK  -  A^AdVjM  -  A.vJk)  W,(A,  -  A,) 

^  (Ad  -  AjvJ  (mV,A2a,  +  KV^Ad  -  MV^AjA^  -  KV^A,] 

in  which  the  terms  involving  matrix  ADd  cancel.  Manipulating  and  simplifying  Eq.  (6.4.11) 
yields 

(a^A^  -  AdAd)vjMVd(Ad  -  A  J  ^  (A^  -  AjvjMVd(AdA2  -  A^Ad)    (6.4.12) 

This  equivalence  is  obviously  satisfied  if  the  eigenvectors  are  mass  orthogonal. 

/// 


(6.4.11) 


PROPERTY  6.4.2  The  updated  FEM,  defined  by  the  original  FEM  and  the  perturbation 
matrices  ADd  and  AKd  computed  from  Eqs.  (6.4.4)  and  (6.4.6),  preserves  the  original  rigid 
body  modes  if  the  measured  eigenvectors  and  the  rigid  body  modes  are  mass  orthogonal. 
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Proof: 

As  discussed  earlier,  a  rigid  body  mode  is  defined  as  a  mode  whose  eigenvalue  is  equal  to 
zero  and  whose  eigenvector  lies  in  the  null  space  of  the  reM  stiffness  matrix.  Hence,  the 
rigid  body  modes  of  the  original  system  are  preserved  in  the  updated  FEM  if  they  lie  in  the 
null  space  of  the  perturbed  stiffness  matrix.  Consider  the  relationship 

e  =  (K  -  AKjvr  (6.4.13) 

where  Vr  is  a  rigid  body  mode  eigenvector.  Clearly,  the  rigid  body  mode  associated  to 
eigenvector  Vr  is  preserved  if  e  =  0.  By  definition,  Vr  is  a  rigid  body  eigenvector  of  the 
original  system,  hence  Eq.  (6.4.13)  can  be  simplified  as 

e  =  -  AK^Vr  (6.4.14) 

Substituting  the  expression  for  AK<j  as  defined  in  Eq.  (6.4.6),  into  Eq.  (6.4.14)  gives 


e  =  -  (BA,  -  BAjH,(BA,  -  BAJ\  (6.4.15) 


Substitution  of  the  expressions  for  B  and  B  into  this  equation  yields 


e  =  -  (BA,  -  BAjHi 


AdAjvjM  +  A,VjK  -  AjAdVjM  -  AjVJk 


Vr  (6.4.16) 


By  using  the  fact  that  Vr  is  a  rigid  body  eigenvector  of  the  original  system  (i.e.  Kvr  =  Q),  Eq. 
(6.3.16)  can  be  simplified  as 

e  =  -  (BA,  -  BAjHk(  A^A^  -  A,Ad  )vJMv,  (6.4.17) 

It  is  clear  from  Eq.  (6.4.17)  that    e    =    0  if  the  rigid  body  mode  Vr  and  the  measured 
eigenvectors  Vj  are  mass  orthogonal  (i.e.  V^Mv^  =  0  ). 


6.4.2     Damage  Extent:  Mass  and  Damping  Properties 

In  this  case  it  is  assumed  that  the  effect  of  the  structural  damage  on  the  stiffness 
properties  is  negligible.  In  this  context,  Eq.  (6.1.1)  is  rewritten  as 
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MVjAj  +  DVjAd  +  KVd  =  AMdV^Aj  +  ADdV^Aj  =  B         (6.4.18) 

By  using  an  approach  similar  to  one  used  in  the  preceding  section,  Eq.  (6.4.18)  and  its 
complex  conjugate  can  be  manipulated  to  yield  the  following  decomposition 

AMdV,(A2A,  -  AdAd)  =  (BA,  -  BAJ  (6.4.19) 

ADdV,(AdAd  -  AdA2)  =  (BAd  -  BA^)  (6.4.20) 

Again  by  applying  the  MRPT  to  the  preceding  equations,  AMd  and  ADd  are  determined  to 
be 

(6.4.21) 


AMj  =  (BAj  -  BAd)Hni(BAd  -  BAJ 


with  Hm  = 


(BAj  -  Ba/v,(a2a,  -  AdAd) 


-1 


AD,  =  ( BAd  -  BAj JHd(BAd  -  BA^) 


with   Hj  = 


T 

(BAd  -  BA2J  V,(AdAd  -  AdA2 


-,-1 


(6.4.22) 


Clearly,  the  perturbation  matrices  AMd  and  ADd  as  defined  by  Eqs.  (6.4.21)  and  (6.4.22)  are 
real. 

PROPERTY  6.4.3  The  perturbation  matrices  AMd  and  ADd,  as  computed  above,  will  be 
symmetric  if  the  measured  eigenvectors,  Vd,  are  stiffness  orthogonal;  i.e.,  the  eigenvectors 
are  orthogonal  with  respect  to  the  original  unperturbed  stiffness  matrix. 

PROPERTY  6.4.4  The  updated  FEM,  defined  by  the  original  FEM  and  the  perturbation 
matrices,  AMd  and  ADd,  preserves  the  original  rigid  body  modes. 

The  proof  of  Property  6.4.4  is  straightforward  since  the  original  stiffness  matrix  is 
unchanged  (see  Property  6.3.2).  The  proof  of  Property  6.4.4  follows  very  much  the  same 
pattern  as  the  proof  of  Property  6.4.2  . 
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6.4.3    Damage  Extent:  Mass  and  Stiffness  Properties 

In  this  problem,  it  is  assumed  that  the  effect  of  the  structural  damage  on  the  damping 
properties  is  negligible.  For  this  situation,  the  general  eigenvalue  problem  defined  in  Eq. 
(6.1.1)  associated  to  this  case  can  be  simplified  as 

MVjA^  +  DVjAd  +  KV,  =  AM,VdA2  +  AK^V,  =  B  (6.4.23) 

Algebraic  manipulafions  of  Eq.  (6.4.23)  and  its  complex  conjugate  yield  the  following 
decomposition 


AM. 


dVd(A3  -  a5)  =  (B  -  B) 
AK,V,(Ad  -  Al)  = 


BAd  -  BA; 

The  perturbation  matrices  AMj  and  AIQi  can  then  be  computed  using  the  MRPT. 


(6.4.24) 
(6.4.25) 


AMd  =  (B  -  B]H„,(B  -  B) 


with  Hm  = 


(B  -  B)\,(a2  -  A^) 


(6.4.26) 


AKh  = 


(BAd  -  BAjjHk^BAd  -  BA^ 


with  Hy.  = 


[baI  -  BAI)\,{aI  -  Aj) 


-,-1 


(6.4.27) 


Note  that  AMj  and  AK^  as  defined  by  Eqs.  (6.4.26)  and  (6.4.27)  are  real. 

PROPERTY  6.4.5  The  perturbation  matrices  AMj  and  AIQ,  as  computed  above,  will  be 
symmetric  if  the  measured  eigenvectors,  Vj,  are  damping  orthogonal;  i.e.,  the  eigenvectors 
are  orthogonal  with  respect  to  the  original  unperturbed  damping  matrix. 


PROPERTY  6.4.6  The  updated  FEM,  defined  by  the  original  FEM  and  the  perturbation 
matrices,  AMj  and  AKd,  preserves  the  original  rigid  body  modes  if  the  measured 
eigenvectors  and  the  rigid  body  modes  are  damping  orthogonal. 
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These  proofs  of  the  above  two  properties  are  not  reported  here.  They  follow  very  much  the 
same  pattern  as  the  proofs  in  Section  6.4.2  . 

6.4.4    Damage  Extent:  Mass.  Damping  and  Stiffness  Properties 

The  eigenvalue  problem  of  a  proportionally  damped  system  with  all  property  matrices 
simultaneously  affected  by  damage  can  be  rearranged  into  the  form 

MV^Aj  +  DVjAd  +  KV,  =  AM,V,A2  +  ADdV^A^  +  AK^Vj  ^  B  (6.4.28) 

The  theory  developed  in  Section  6.3.3  can  be  expanded  to  address  this  particular  problem. 
The  cross-orthogonality  relationships  associated  with  this  type  of  structures  are 

Vj(M  -  AMJV,  =  Ip,p  (6.4.29a) 

Vj(D  -  ADjVd  =  diag(2Cd  Wd,,    •  •  •  ,2^^,^^^)  =  ^d  (6-4.2%) 

Vj(K  -  AKjVj  =  diag(a)d_2,    .  .  .  ,o),^2J  =  q^  ^^^^Qc) 

Notice  that  the  cross-orthogonality  relationships  in  Eqs.  (6.4.29a)  and  (6.4.29c)  are  exactly 
the  same  as  the  ones  associated  with  undamped  systems  reported  in  Eqs.  (6.3.16a)  and 
(6.3. 16b).  As  before,  these  cross-orthogonality  conditions  can  also  be  rearranged  as 

VjAMjVd    =    VjMVd    -    Ipxp    =    VjB^  (6.4.30a) 

VjAD,  Vd    =    VjDV,    -    Sd    =    VjBd  (6.4.30b) 

VjAK^Vd    =    VjKV,    -    Q,    =    VjBk  (6.4.30c) 

Following  the  exact  same  argument  discussed  for  undamped  systems  in  Section  6.3.3,  an  nxp 
matrix  P  that  satisfies  the  relation, 

p(vJb)     =    B  (6.4.31) 

is  sought,  where  B  is  computed  using  Eq.  (6.4. 18)  and  (vJb  j  is  a  pxp  matrix.  Although  B 
is  a  complex  matrix,  the  nxp  matrix  P  is  real,  since  Vj  is  real.  Hence,  for  computational 
efficiency,  matrix  P  can  computed  from 
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P    =    Br(vTBr)  (6.4.32) 

where  Br  is  the  real  part  of  B.  In  Eq.  (6.4.32),  it  is  assumed  that  matrix  (vJb  j  is  invertible. 
With  P  computed,  the  next  step  is  to  determine  the  decomposed  damage  vectors  that  indicate 
the  effects  of  the  damage  on  the  mass,  damping  and  stiffness  matrices, 

Bm    =    AMVd    =    P(VjMVd    -    Ipxp)  (6.4.33a) 

B,    ^    ADV.    =    P(VTdV.    -    2.)  (6.4.33b) 


i^    =    AKV,    =    P(Vj 


B.    =    AKV.    =    P(VTkV.    -    Q.)  (6.4.33c) 


The  minimum  rank  perturbation  theory  (MRPT),  as  formulated  in  Proposition  6. 1 ,  can  again 
be  appHed  to  determine  the  perturbation  matrices,  AMj,  ADj  and  AIQi,  as 

AMj  =  Bm  (B^Vd)-'  B;^  (6.4.34a) 

ADd  =  B,  (BjVj)-'  Bj  (6.4.34b) 

AK^  =  Bk  (B^V^)-'  Bj  (6.4.34c) 

Note  that  the  matrices  B^Vj  ,  bJv^  and  bJv^  are  pxp  matrices  that  are  invertible  if  Bm, 
Bd  and  B^  are  of  full  rank.  As  in  all  other  cases  already  studied.  Proposition  6.2  can  be  used 
to  deal  with  the  situation  when  any  one  of  these  matrices  are  rank  deficient.  The  cumulative 
damage  location  vector  associated  to  B^  and  Bk,  defined  in  Eqs.  6.3.21,  are  also  applicable 
to  this  problem.  An  additional  cumulative  damage  vector  associated  to  the  perturbations  in 
the  damping  properties  can  be  similarly  defined  as 

where  d^  is  the  i'*'  column  of  matrix  Bj. 
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PROPERTY  6.4.5  The  perturbation  matrices  (AMd,  ADj,  AK^)  as  computed  above  will  be 
symmetric. 

Proof: 

This  is  true  since  these  perturbation  matrices  are  constrained  to  satisfy  the  relationships  in 

Eqs.  (6.4.30).  Notice  that  the  right  hand  sides  of  Eqs.  (6.4.30)  are  symmetric. 


6.5    Damage  Extent:  Nonproportionally  Damped  Structures 

As  reported  earlier  in  the  beginning  of  the  chapter,  the  eigenvalue  problem  associated 
with  the  dynamic  differential  equation  of  motion  of  a  damaged  nonproportionally  damped 
system  can  be  written  as 

MV,A2  +  DVjA,  +  KV,  =  AM^VjA^  +  AD.V.A^  +  AK,V,  ^  B  (6.5.1) 

where  the  matrices  M,  D,  K,  AMj,  ADj,  AK^,  Vj,  and  A^j  have  the  same  definitions  as  in 
the  previous  sections.  A  number  of  cross-orthogonality  relationships  associated  to  general 
nonproportionally  damped  systems  have  been  reported  in  Chapter  2.  The 
cross-orthogonality  conditions  that  are  relevant  to  the  following  formulations  are 


VdAd 


VdAd 


VdAd 


'(M  -  AMj)  0 

0  -  (K  -  AK,) 

■(D  -  AD,)  (K  _  AKj' 

(K  -  AKj  0 

0  (m  -  amJ' 

(M  -  AMd)  (D  -  AD,) 


irvdAd" 

"o" 

J[  V,  _ 

0 

[VdAdj 

r 

0' 

[  Vd 

■ 

0 

IrvdAd- 

"o" 

Vd 

0 

J 

(6.5.2) 


(6.5.3) 


(6.5.4) 


In  the  above  equations,  the  notation  [  J  denotes  the  complex  conjugate  transpose  operator. 
Expansion  of  the  above  three  equations  yield,  respectively,  the  following  relations, 

AdV^(M  -  AM,)V,A,  -  V^(K  -  AK,)V,  =  [0]  (6.5.5) 
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AdVd(D  -  ADjVjAd  +  AdVd(K  -  AKJv,  +  Vd(K  -  AK^jV^Ad  =  [0]    (6.5.6) 

Vd(M  -  AM)VdAd  +  AdVd(M  -  AM^V^  +  Vd(D  -  AD)^^  =  [0]    (6.5.7) 

These  equations  can  be  further  manipulated  to  become,  respectively, 

W,  ^  AdVd  M  VdAd  -  VS  K  Vj 

(6.5.8) 

=  A^Vd  AM,  W,A,  -  Y*d  AKd  V, 

W2  =  AdVd  D  VjAd  +  AdVd  K  Vj  +  Vd  K  VjA, 

(6.5.9) 

=  A^V^  AD,  V,A,  +  AdVd  AK,  V,  +  V^  AK,  V,A, 

W3  ^  Vd  M  V,A,  +  AdVd  M  V,  +  Vd  D  V, 

(6.5.10) 
=  V^  AM,  V,A,  +  A^V^  AM,  V,  +  V^  AD,  V, 

Notice  that  Wj ,  W2  and  W3  can  be  computed  using  the  equations  in  terms  of  known  matrices 
M,  D,  K,  Vd,  and  Ad. 

The  above  cross-orthogonality  conditions  are  derived  from  two  different  state  space 
representations  of  the  equation  of  motion  shown  in  Eq.  (6.5. 1)  (see  Chapter  2).  From  each 
representation,  two  cross-orthogonality  condition  were  derived.  The  two  sets  provide 
exactly  the  same  information  about  Eq.  (6.5. 1)  but  in  a  different  format.  Hence,  only  two  of 
the  three  cross-orthogonahty  relations  are  linearly  independent.  In  fact,  Eq.  (6.5.8)  can  be 
generated  by  subtracting  Eq.  (6.5.9)  from  (6.5.10).  Any  set  of  two  linearly  independent 
cross-orthogonality  relationships  can  be  viewed  as  a  decoupling  of  the  equation  of  motion 
(Eq.  (6.5.1)).  Hence,  for  the  present  problem,  there  are  only  two  distinct  (linearly 
independent)  matrix  equations  for  the  three  unknown  perturbation  matrices  (AMj,  ADj, 
AKd).  Naturally,  if  the  problem  consists  of  solving  for  all  three  unknown  matrices,  there  will 
exist  an  infinite  number  of  solutions.  Additional  constraints  are  hence  needed  to  get  a  better 
defined  problem.  In  this  formulation,  it  is  assumed  that  the  damage  only  affects  two  of  the 
three  property  matrices.  Mathematically,  this  can  be  translated  as  a  constraint  on  one  of  three 
perturbation  matrices  (either  AMd  or,  ADd  or,  AIQ)  to  be  equal  to  zero. 
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6.5.1     Damage  Extent:  Damping  and  Stiffness  Properties 

In  this  case,  it  is  assumed  that  the  effects  of  the  damage  on  the  mass  properties  are 
negligible  (AMj  =  0).  The  corresponding  eigenvalue  equation  and  the  cross-orthogonality 
conditions  relevant  to  this  problem  are 

MV^Aj  +  DVjAd  +  KV,  =  ADjVdAd  +  AK^V,  =  B  (6.5.11) 

Vd  AKj  Vd  =  -  Wi  (6.5.12) 

Vd  ADd  Vj  =  W3  (6.5.13) 

The  approach  used  in  this  formulation  follows  very  much  the  same  pattern  as  the  one 
developed  in  Section  3.3.3.  The  basic  ideals  to  generate  fromEqs.  (6.5.1 1-13)  problems  that 
can  take  the  advantage  of  the  already  existing  MRPT  concept.  Based  on  the  same  argument 
as  in  Section  3.3.3,  a  matrix  P  G  C""?  is  sought  such  that 

p(VdB)    =    B  (6.5.14) 

It  is  Assumed  that  (vJb)  is  invertible;  then  matrix  P  is  given  by 

P    =    B  (VdB)  (6.5.15) 
Matrix  P  can  then  be  used  in  Eqs.  (6.5.11)  and  (6.5.12)  to  yield 

AKd  Vd  =  -  PW,  (6.5.16a) 

ADd  Vj  =  PW3  (6.5.16b) 

In  order  to  force  the  perturbation  matrices  ADd  and  AKd  to  be  real,  the  complex  conjugates 
of  Eqs.  (6.5.16a)  and  (6.5.16b)  are  introduced. 


AKj  Vd  =  -  PWi  (6.5.17a) 


ADd  ^d  =  PW3  (6.5.17b) 

A  combination  of  Eqs.  (6.5.16a,  b)  and  (6.5. 17a,  b)  implies 
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AK 


d 

AD 


(VdJR   I   (Vd) 


(PW,)r   I   (PW,) 

(PW3)  I  (PW3) 


(6.5.18a) 
(6.5.18b) 


where  the  notation  (X)^  and  (X)j  indicate,  respectively,  the  real  and  imaginary  part  of  matrix 
X.  The  MRPT  can  now  be  applied  to  determine  the  perturbation  matrices  ADj  and  AKd. 


AD,  =  F,  (fJuJ"    f; 


AKd  =  F,  (F^U,)"    F' 


(6.5.19a) 
(6.5.19b) 


where 


Ud  = 


(VcIr   '   (Vd),] 
(PW,)r   I   (PW,),] 
(PW3)^    I    (PW3) 


Fd    = 


Remember  that  Wi  and  W3  can  be  computed  in  terms  of  known  matrices  as  shown  in  Eqs. 
(6.5.8)  and  (6.5.10).  The  computations  in  Eqs.  (6.5.19a)  and  (6.5.19b)  require  that  (fJU,) 

and  (fJUj j  be  invertible.  As  discussed  in  the  earlier  sections.  Proposition  6.2  can  be  used 
when  these  inverse  requirements  are  not  meet. 

6.5.2     Damage  Extent:  Mass  and  Damping  Properties 

Here,  it  is  assumed  that  the  damage  effects  on  the  stiffness  properties  are  negligible 
(AIQi  =  0).  By  considering  this  assumption,  the  eigenvalue  equation  associated  with  this 
problem  is  found  to  be 

MV,A2  +  DVjAd  +  KVd  =  AM^V^Aj  +  AD^V^Ad  =  B         (6.5.20) 

The  cross-orthogonality  conditions  relevant  to  the  present  development  are  obtained  from 
Eqs.  (6.5.8)  and  (6.5.9)  by  setting  AIQ  to  zero. 


AdVd  AM,  VjAj  =  Wi 


(6.5.21) 
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AdVd  AD,  V,A,  =  W2 


(6.5.22) 


The  solution  procedure  to  the  problem  of  determining  AMj  and  ADj  is  exactly  the  same  as 
the  one  formulated  in  the  preceding  section.  First,  equations  (6.5.21)  and  (6.5.22)  are  cast 
in  the  following  form. 


AM,  V,A,  =  PW, 


AD,  V,A,  =  PW2 


where      P    =    B  (AdVdB) 


(6.5.23a) 
(6.5.23b) 


-1 


The  above  two  equations  involve  complex  matrices.  To  guarantee  that  the  perturbation 
matrices  AMd  and  ADd  are  real,  Eqs.  (6.5.23a)  and  (6.5.23b)  are  used  in  conjunction  with 
their  corresponding  complex  conjugates.  The  combinations  of  the  two  sets  of  complex 
conjugate  pairs  yield 


AM. 


AD. 


(VdAd)R   I   (VdA,), 
'(VdAd)R   I   (VdAd), 


(PW,)r   I   (PW,)j]  (6.5.24a) 

"(PWj   I   (PWjJ  (6.5.24b) 


Finally,  the  perturbation  matrices  AMd  and  ADd  ^e,  respectively,  determined  from  Eqs. 
(6.5.24a)  and  (6.5.24b)  by  using  the  MRPT. 


-1 


AM,  =  Fm  (F^U,)       Fj 
ADd  =  F,  (fJu,)"'fJ 


(6.5.25a) 
(6.5.25b) 


where 


Ud  = 


Fh 


(VdAdJR   I   (VdAd)J 
(PWi)r   I   (PW,)j 
(PW2)r   I   (PW2),] 


p.-.:: 
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Again,  the  computations  in  Eq.  (6.5.25a)  and  (6.5.25b)  are  feasible  if  and  only  if  (PmUj)  and 
(fJUj)  are  both  invertible. 

6.5.3    Damage  Extent:  Mass  and  Stiffness  Properties 

In  some  situations,  structural  damage  is  such  that  it  only  affects  the  mass  and  stiffness 
properties  of  the  system.  The  following  expression  are  the  results  of  setting  the  damping 
perturbation  to  zero  (ADj  =  0)  in  Eqs.  (6.5.1),  (6.5.9)  and  (6.5.10). 

MVjA^  +  DVjAd  +  KVd  =  AM^VjA^  +  AK^Vj  ^  B  (6.5.26) 

AdVd  AKd  V,  +  Vd  AKd  VjA,  =  W2  (6.5.27) 

AdVd  AMj  Vd  +  Vd  AMj  VjAd  =  W3  (6.5.28) 

Equation  (6.5.26)  is  the  rearranged  eigenvalue  problem  associated  with  the  present  case. 
The  expressions  in  Eqs.  (6.5.12)  and  (6.5.13)  are  the  "rearranged"  cross-orthogonality 
conditions  that  are  relevant  to  the  forthcoming  formulation.  First,  Eqs.  (6.5.27)  and  (6.5.28) 
are  conveniently  simplified  to  the  standard  format  encountered  in  the  previous  two  sections. 
Since  Ad  and  Ad  are  diagonal  matrices,  Eqs.  (6.5.27)  can  be  written  as 

Vd  AKd  Vd   =  Q2  (6.5.29) 

where      q2_^   =  v^2,  I  (^1,  +  ^dj 
The  variables  q2  and  Wj  are  the  components  in  the  i'*'  row  and  j*  column  of  matrix  Qk  and 
W2,  respectively;  "k^  is  the  i'^  diagonal  component  of  the  eigenvalue  matrix  Aj.  Similarly, 
Eq.  (6.5.28)  can  also  be  written  as 

Vd  AMd  Vj   =  Q3  (6.5.30) 

where      qg,^    =   W3,.  /  \\X  +  XdJ 
Equations  (6.5.29)  and  (6.5.30)  can  then  be  transformed  as 


F- 
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AK,  V,  =  PQ^ 


AM,  Vd  =  PQ3 


(6.5.31a) 


(6.5.31b) 


where      P    =    B  (VdB) 

Notice  that  the  above  computations  require  that  matrix  (VdB)  be  invertible.  The  sought 

perturbation  matrices,  AMj  and  AK^,  are  ensured  to  be  real  by  supplementing  Eqs.  (6.5.3  la) 
and  (6.5.31b)  with  their  complex  conjugates.  The  complex  conjugate  pairs  are  then 
combined  together  in  the  following  fashion 


AKd 
AMh 


(PQ2)r    I    (PQ2),] 
(PQ3)    I    (PQ3),] 


(6.5.32a) 
(6.5.32b) 


At  this  point,  the  MRPT  is  applied  to  Eqs.  (6.5.32a)  and  (6.5.32b)  and  the  perturbation 
matrices,  AMj  and  AK^i,  are  determined 

-1 


AKd  =  F,  (FjuJ       F' 


AM,  =  Fn,  (pTuJ       F 


(6.5.33a) 


(6.5.33b) 


where 


Ud  = 


Fk    = 


Fn,      = 


(VdJR      I      (Vd) 

(PQ2)r    I    (PQ2)i] 

(pQ3)r  I  (pQb),] 


Again,  it  is  assumed  that  Proposition  6.2  has  been  applied,  if  necessary,  to  ensure  that  the 
inverses  of  the  matrices  (FmU,)  and  (fJu,)  exist. 

PROPERTY  6.5.1  In  the  above  three  damage  configurations  considered  for 
nonproportionally  damped  structures,  the  computed  perturbation  matrices  (AMj,  ADj  or 
AKd  )  from  the  MRPT  will  be  symmetric. 


*      f 
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Proof: 

This  is  true  since  tlie  developed  formulations  require  the  perturbation  matrices  along 
with  the  "healthy"  FEM  to  satisfy  the  cross  orthogonality  conditions  in  Eqs.  (6.5.2-4).  These 
cross-orthogonality  conditions  are  satisfied  if  and  only  if  the  corresponding  state  matrices 
are  symmetric;  that  is 

-iT 

(6.5.34) 


(M  -  AM^)  0 

0  -  (K  -  AK,) 


(M  -  AMj)  0 

0  -  (K  -  AKj 


■(D  -  AD,)  (K  -  AK,j 

(K  -  AKj  0 

0  (m  -  amJ" 

(M  -  AM,)  (D  -  ADj 


■(D  -  AD,)  (K  -  AK,)" 
(K  -  AK,)  0 

0  (M  -  AM,)' 
(M  -  AM,)       (D  -  AD,) 


(6.5.35) 


(6.5.36) 


Since  the  matrices  M,  D,  K  are  known  to  be  symmetric,  the  above  relationships  are  satisfied 
if  and  only  if  the  perturbation  matrices  AMj,  ADj  and  AK<j  are  symmetric. 


6.6    Practical  Issues 


In  all  cases  discussed  in  the  previous  sections  of  this  Chapter,  the  damage  extent  problem 
boiled  down  to  solving  problems  of  the  form 

AX    =    Y  (6.6.1) 

in  which  A  is  an  unknown  nxn  matrix  representing  one  of  the  perturbation  matrices  (AMj, 
ADd  and  AK^)  and  X  and  Y  are  known  nxp  complex  matrices.  Matrix  X  is  a  function  of  the 
damaged  structure  eigenvector  (Vj)  and  eigenvalue  (Ad)  matrices.  Recall  that  p  represents 
the  number  of  modes  and  n  is  the  total  number  of  DOFs  in  the  structure's  FEM.  As  discussed 
in  Chapter  5,  the  columns  of  matrix  Y  provide  useful  information  in  determining  the  location 
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of  the  structural  damage.  For  each  problem,  a  cumulative  damage  vector  associated  to  the 
columns  of  matrix  Y  can  be  defined  as 


1  p    lyJ 


Ydi 


Recall  that  the  DOFs  of  the  FEM  affected  by  the  damage  are  associated  with  the  components 
of  vector  d  which  are  substantially  larger  than  the  other  components  of  d. 

Considering  the  present  generic  problem,  the  extent  problem  consists  of  determining  a 
symmetric  real  matrix  A  that  satisfies  Eq.  (6.6.1).  A  solution  to  this  problem  based  on  the 
Minimum  Rank  Perturbation  Theory  (MRPT)  was  developed  in  Proposition  6.1.  The 
solution  technique  presented  in  Proposition  6.1  requires  the  pxp  matrix  (Y^X)  to  be 

invertible.  Proposition  6.2  was  then  presented  to  deal  with  the  case  when  matrix  (Y^X)  is  not 
invertible.  In  the  forthcoming  discussions,  two  strategies  to  improve  the  extent  calculation 
are  presented. 

6.6.1     The  Concept  of  "Best"  Modes 

Note  that  for  the  sake  of  generality,  the  forthcoming  discussion  is  based  on  the  generic 
problem  shown  in  Eq.  (6.6.1);  however  it  should  be  viewed  as  representative  of  all  extent 
problems  discussed  in  the  earlier  sections  of  this  chapter. 

Although  Proposition  6.2  deals  with  noise  free  measurements,  its  practical  usefulness  is 
for  the  case  when  the  eigenvalue  /  eigenvector  measurements  are  corrupted  by  noise  (either 
through  actual  measurement  error  or  eigenvector  expansion  induced  errors).  In  general, 
matrix  Y  will  typically  be  of  full  column  rank  due  to  the  presence  of  noise  ( rank(Y)=p ),  even 
though  the  actual  damage  has  only  induced  a  rank  r  (  r  <  p  )  change  to  the  FEM  property 
matrix  (M,  D  or  K)  associated  to  A  ( AMj,  ADd  or  AK^).  Application  of  the  extent  algorithm 
at  this  point  will  result  in  a  rank  r  perturbation  to  the  property  matrix  A.  A  better  estimate  of 
the  true  extent  of  the  damage  can  be  obtained  by  estimating  the  "true"  rank  of  Y  (and  thus  A), 
where   the   true   rank   of  Y   is   defined   by   the   limiting   case   when   there   are   no 
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eigenvalue/eigenvector  measurement  errors.  There  are  two  possible  techniques  for 
estimating  the  "true"  rank  of  Y.  The  first  technique,  which  amounts  to  a  brute-force 
numerical  computation,  is  to  calculate  the  Singular  Value  Decomposition  (SVD)  of  Y. 
Various  rank  criteria  measures,  as  formulated  by  Juang  and  Pappa  (1985)  in  their  work  on 
system  realization  theory,  in  conjunction  with  engineering  judgement,  can  be  used  to  arrive 
at  an  estimate  of  the  true  rank  of  Y  and,  hence,  A.  However,  this  would  greatly  increase  the 
computational  burden  of  the  extent  algorithm,  and  may  not  even  be  practical  if  the  number  of 
DOFs  of  the  original  FEM  is  large.  A  computationally  inexpensive  technique  for  estimating 
the  true  rank  of  Y  is  to  make  use  of  the  damage  location  results.  The  cumulative  damage 
vector  associated  to  the  perturbation  matrix  A  defined  in  Eq.  (6.6.2)  indicates  which 
structural  DOFs  have  been  directly  affected  by  the  damage.  Using  connectivity  information 
from  the  original  FEM  in  conjunction  with  knowledge  of  the  damaged  DOFs  allows  for  the 
determination  of  which  elements  of  the  FEM  have  been  damaged.  The  "true"  rank  of  Y  can 
be  estimated  by  adding  up  the  ranks  of  the  element  matrices  associated  with  each  damaged 
element. 

Once  the  true  rank  of  Y  has  been  estimated,  corresponding  columns  of  X  and  Y  are 
eliminated  such  that  the  actual  rank  of  B  is  set  to  the  estimated  true  rank  of  Y  In  Proposition 
6.2,  it  was  shown  that  with  noise  free  measurements  it  does  not  matter  which  r  modes  from 
the  set  of  p  measured  modes  are  used  in  the  extent  calculation.  However,  when  the 
measurements  are  corrupted  by  noise,  it  is  best  to  use  those  modes  which  best  reflect  the 
damage  state  as  indicated  by  the  cumulative  damage  vector  Eq.  (6.6.2).  As  will  be  illustrated 
in  the  example  problems  presented  in  Chapter  7,  the  concept  discussed  here  improves  the 
damage  extent  estimation.  This  is  especially  true  when  the  rank  of  the  perturbation  matrix  (r) 
due  to  the  damage  is  less  than  the  number  of  measured  modes  (p).  The  reason  is  twofold. 
First,  the  rank  of  the  computed  perturbation  matrix  is  consistent  with  the  rank  of  the  actual 
perturbation  matrix.  Second,  the  measurement  errors  in  the  experimental  modes  that  do  not 
reflect  the  actual  state  of  the  damage  are  not  included  in  the  computation. 
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6.6.2     Application  of  the  Eigenvector  Filtering  Algorithm 

Another  technique  that  can  be  used  to  improved  the  damage  extent  estimation  is  the 
eigenvector  filtering  algorithm  developed  in  Chapter  5.  The  modes  that  should  be  filtered 
and  used  in  the  extent  computation  are  the  ones  that  reflect  the  overall  structure  cumulative 
damage  vector  defined  in  Eqs.  (5.4.1)  or  (5.4.2)  of  Chapter  5.  The  characteristics  of  the 
eigenvector  filtering  algorithm  are  illustrated  in  Chapter  7.  It  should  be  pointed  out  that  it  is 
not  a  requirement  of  the  extent  algorithm  to  use  the  eigenvector  filtering  algorithm,  but 
certainly  it  does  allow  for  more  engineering  judgement  to  enter  into  the  extent  calculation. 
The  choice  of  whether  or  not  to  use  the  filtering  algorithm  is  dependent  on  (i)  the  errors  in  the 
eigenvectors  as  seen  from  the  damage  vectors,  (ii)  the  requirement  on  accuracy  of  the  extent 
calculation  and  (iii)  the  size  of  the  FEM,  which  defines  the  additional  computational  burden. 

6.7    Summary 

A  computationally  attractive  algorithm  was  developed  to  provide  an  insight  to  the  extent 
of  structural  damage.  The  developed  theory,  called  the  Minimum  Rank  Perturbation  Theory 
(MRPT),  makes  use  of  an  existing  finite  element  model  of  the  "healthy"  structure  and  a 
subset  of  experimentally  measured  modal  properties  of  the  "damaged"  structure.  The  MRPT 
formulation  is  consistent  with  the  effects  of  many  classes  of  structural  damage  on  a  finite 
element  model.  The  MRPT  was  applied  to  assess  the  damage  in  undamped,  proportionally 
damped  and  nonproportionally  damped  models.  For  each  type  of  model,  several  different 
damage  scenarios  were  considered.  The  performance  of  the  MRPT  will  be  investigated  in 
Chapter  7  along  with  the  subspace  rotation  damage  location  algorithm  using  both  simulated 
and  actual  experimental  data.  As  will  be  shown  in  Chapter  7,  the  MRPT  can  also  be  used  to 
refine  finite  element  models. 


CHAPTER  7 

VALIDATION  AND  ASSESSMENT  OF  THE  SUBSPACE  ROTATION  ALGORITHM 

AND  THE  MINIMUM  RANK  PERTURBATION  THEORY 


7.1    Introduction 

In  this  chapter,  the  characteristics  of  the  subspace  rotation  algorithm  (Chapter  5)  and  the 
minimum  rank  perturbation  theory  (Chapter  6)  are  illustrated  by  using  both  simulated  and 
actual  experimental  data.  The  main  objective  is  to  show  that  the  algorithms  are  suited  to 
assess  damage  in  structures.  Additionally,  it  is  shown  that  the  minimum  rank  perturbation 
theory  formulated  in  Chapter  6  is  also  suited  to  update  finite  element  models.  The  computer 
simulated  example  problems  provide  a  controlled  setting  in  which  key  points  made 
throughout  the  formulation  of  both  algorithms  are  emphasized.  Actual  experimental 
examples  are  used  to  demonstrate  the  practicality  of  the  algorithms  in  handling  "real-life" 
systems. 

7.2    Kabe's  Problem 

7.2.1     Problem  Description 

Kabe's  eight  degree  of  freedom  spring-mass  problem  is  shown  in  Figure  7.1,  which 
includes  the  stiffness  and  mass  values  for  the  exact  model.  The  main  objective  of  this 
problem  is  to  illustrate  the  damage  location  algorithm  formulated  in  Chapter  5.  The 
Subspace  Rotation  algorithm  direct  method  is  compared  to  the  Angle  Perturbation  Method. 
Key  points  discussed  in  the  formulation  will  be  highlighted.   Application  of  the  damage 
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extent  algorithm  of  Chapter  6  is  also  addressed.    Both  algorithms  will  be  evaluated  using 
noisy  modal  measurements  in  an  attempt  to  reproduce  a  practical  testing  situation. 
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mi  =  0.001      m8  =  0.002     mj  =  1.0     j  =  2, ...  ,7 
ki  =  1000     k2=10     k3  =  900     k4=100     k5=1.5     k6  =  2.0 

Figure  7. 1    Kabe's  Problem. 

A  variation  of  Kabe's  original  problem  is  used  here.  Rather  than  the  standard  model 
commonly  used,  which  has  incorrect  values  for  ail  of  the  connecting  springs,  only  a  single 
spring  constant  is  changed.  This  is  reflective  of  the  fact  that  damage  may  occur  as  a  large 
local  change  in  the  stiffness  of  a  structural  member.  In  this  problem,  the  initial  model  is 
incorrect  for  the  spring  between  masses  seven  (7)  and  eight  (8).  Damage  will  be  modeled 
as  a  change  in  spring  constant  from  1.5  to  0.1.  Changing  the  spring  in  this  fashion  causes 
the  damage  to  occur  at  degrees  of  freedom  which  have  a  small  ||  Z^  ||  in  comparison  to  other 
DOFs.  Thus,  it  would  be  expected  that  damage  would  best  be  located  using  the  Angle 
Perturbation  method  as  defined  in  Eq.  (5.3. 1).  In  this  problem,  it  is  assumed  that  only  the 
first  two  modes  of  vibration  are  measured.  Note  that  the  damage  does  not  affect  the  mass 
properties  of  the  structure. 
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7.2.2     Damage  Location 

Figures  7.2  and  7.3  show,  respectively,  the  damage  location  results  of  using  the 
Subspace  Rotation  direct  (dj)  method  and  the  Angle  Perturbation  Method  (ctj).  Figure  7.4 
provides  the  damage  vector  as  defined  by  Lin's  algorithm  (1990).  These  results  were 
generated  using  only  the  fundamental  mode  of  vibration  eigendata.  In  these  figures,  the 
x-coordinate  represents  the  elements  of  the  damage  vector  di  (or  aj),  where  the  i^*'  element 
of  dj  (aj)  corresponds  to  the  i'^  degree  of  freedom. 

The  upper  left  plot  represents  the  damage  vector  if  the  exact  perturbation  matrices  are 
used.  The  upper  right  plot  corresponds  to  the  case  where  the  exact  eigenvector  information 
is  provided  to  the  algorithm.  The  lower  left  and  lower  right  plots  correspond  to  the  case 
where  the  exact  eigenvector  has  been  corrupted  with  five  and  ten  percent  random  noise, 
respectively.  In  practice,  the  noise  in  eigenvector  information  could  be  due  to  both 
measurement  and/or  expansion  induced  errors. 

As  shown  in  Figures  (7.2)  and  (7.3),  the  location  algorithm  is  able  to  exactly  locate  the 
correct  damage  when  provided  with  only  a  single  correct  eigenvalue/eigenvector  pair  as 
guaranteed  by  Eq.  (5.2.4).  From  the  lower  plots,  it  is  seen  that  the  algorithm  experiences 
a  degradation  in  damage  detection  capability  as  the  error  in  the  eigenvector  is  increased.  In 
fact,  when  presented  with  10%  eigenvector  noise,  the  direct  subspace  rotation  method 
incorrectly  identifies  the  damage.  However,  when  using  the  angle  perturbation  method  the 
location  algorithm  is  able  to  discern  damage  with  a  10%  noise  level.  As  shown  in  Figure 
7.4,  Lin's  algorithm  is  unable  to  discern  damage  in  this  particular  test  problem. 

Figures  7.5  and  7.6  show  the  results  of  using  the  angle  perturbation  method  and  Lin's 
algorithm  when  provided  with  complete  first  and  second  mode  eigenvalue/eigenvector 
information.  Again,  the  angle  perturbation  method  is  able  to  cleanly  detect  damage  even 
when  presented  with  eigenvector  information  corrupted  with  ten  percent  noise.  Again,  Lin's 
algorithm  is  unable  to  discern  damage. 
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Figure  7.2     Kabe's  Problem:  Damage  Location  Results  using  the  Subspace 
Rotation  Direct  Method  with  the  Eigendata  of  the  1st  Mode. 
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Figure  7.3     Kabe's  Problem:  Damage  Location  Results  using  the  Angle 
Perturbation  Method  with  the  Eigendata  of  the  1st  Mode. 
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Figure  7.4     Kabe's  Problem:  Damage  Location  Results  using  Lin's 
Algorithm  with  the  Eigendata  of  the  1st  Mode. 
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Figure  7.5     Kabe's  Problem:  Damage  Location  Results  using  the  Angle  Perturbation 
Method  with  the  Eigendata  of  the  1st  and  2nd  Modes. 
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Figure  7.6     Kabe's  Problem:  Damage  Location  Results  using  Lin's 
Algorithm  with  the  Eigendata  of  the  1st  and  2nd  Modes. 


7.2.3     Damage  Extent 


In  the  previous  two-mode  noisy  test  case,  inspection  of  the  cumulative  damage  vector 
clearly  indicates  that  the  7th  and  8th  degrees  of  freedom  have  experienced  damage. 
Inspecting  the  structural  connectivity,  it  can  be  deduced  that  the  spring  between  mass  7  and 
8  has  been  damaged.  The  rank  of  the  single  spring  "element"  matrix  is  one.  Thus,  the  rank 
of  the  "true"  perturbation  to  the  stiffness  matrix  due  to  damage  is  one.  From  Proposition  6. 1 
and  6.2,  it  is  clear  that  only  experimental  data  from  one  mode  of  vibration  should  be  used 
to  compute  the  extent  of  damage.  In  the  noisy  cases,  as  discussed  in  Section  6.6. 1 ,  the  mode 
that  should  be  used  in  the  extent  calculation  is  the  one  associated  with  a  dj  that  most  cleanly 
demonstrates  the  damage  location  shown  in  Figure  7.5.  A  simple  inspection  of  the  individual 
dj  of  the  two  measured  modes  indicates  that  mode  2  provides  the  best  insight  into  the  state 
of  damage.  Figure  7.7  presents  element-by-element  stiffness  perturbation  matrix  results 
from  the  application  of  the  extent  algorithm,  formulated  in  Section  6.3.3,  using  mode  2 
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eigendata.  The  x-coordinates  on  each  plot  are  the  indices  of  a  column  vector  constructed  by 
storing  the  upper  triangular  portion  of  the  perturbation  stiffness  matrix  (AFQj)  in  a  column 
vector.  The  y-coordinate  on  all  plots  consists  of  the  difference  between  the  updated  stiffness 
matrix  and  the  original  stiffness  matrix  (i.e.  AK<i).  In  the  noise  free  case,  it  should  be  noted 
that  the  exact  damage  would  be  computed  using  any  single  set  of  eigenvalue/eigenvector 
measurement.  The  extent  algorithm  experiences  a  modest  degradation  as  the  measurements 
are  corrupted  by  noise. 
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Figure  7.7     Kabe's  Problem:  Damage  Extent  Results  using  the  MRPT 
with  the  Eigendata  of  Mode  2. 

To  contrast  the  proposed  minimum  rank  perturbation  approach  to  the  optimal  matrix 
update  approaches,  consider  the  application  to  the  current  problem  of  the  commonly  used 
method  formulated  by  Baruch  and  Bar  Itzhack  (1978).  Although  this  algorithm  was 
developed  from  a  model  refinement  viewpoint,  it  has  been  investigated  for  possible  uses  in 
damage  detection  (Smith,  1992).  It  is  evident  from  Figure  7.8  that  this  algorithm  is  unable 
to  discern  the  damage  for  this  particular  example  when  provided  with  (i)  the  second  mode's 
exact  eigendata,  or  (ii)  the  first  six  modes'  exact  eigendata.  In  fact,  it  can  be  seen  in  the  upper 
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Figure  7.8     Kabe's  Problem:  Damage  Extent  Results  using  Baruch's  Method. 

right  plot  of  Figure  7.8  that,  consistent  with  the  minimum  Frobenius  norm  change 
formulation,  the  approach  tends  to  "smear"  the  changes  throughout  the  entire  stiffness 
matrix  instead  of  localizing  them  at  the  damage  location.  In  this  problem,  this  method  is  only 
able  to  cleanly  ascertain  the  damage  when  provided  with  the  exact  eigendata  of  the  first 
seven  modes. 

7.3    Damage  Detection:  Fifty-Bay  Two-Dimensional  Truss:  Undamped  FEM 


7.3.1     Problem  Description 

The  50-bay,  2-dimensional  truss  used  in  this  example  is  shown  in  Figure  7.9.  The 
geometric  and  material  properties  of  the  truss  are  also  given  in  the  figure.  Each  truss  member 
was  modeled  as  a  rod  element.  The  finite  element  model  of  the  structure  has  201  degrees 
of  freedom  and  is  undamped.  In  this  example,  damage  is  simulated  by  reducing  the  Young's 
modulus  of  two  members.  One  of  the  damaged  members  is  the  upper  longeron  of  the  third 
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bay.  In  this  member,  the  Young's  modulus  is  decreased  from  E=2.0xl0^^  Pa  (29x10^  psi) 
to  6.895x10^  Pa  (1x10^  psi).  The  other  member  is  the  lower  longeron  of  bay  forty.  This 
member  is  subjected  to  a  complete  loss  of  stiffness  (Young's  modulus  equal  to  zero).  It  is 
assumed  that  only  the  first  ten  modes  of  vibration  are  measured. 


Bay  1 


// 


// 


Bay  50 


'  ■ 


E  =  2.0x1011  N/m2  (29x10^  psi)  p  =  7.833xl03  kg/m^  (0.283  Ib/in^) 

A  =  6.452X  10-^  m2  (1.0  in^)  L  =  0.1270  m  (5.0  in) 

Figure  7.9     Fifty-Bay  Two-Dimensional  Truss. 

7.3.2     Damage  Location 

The  first  task  in  this  analysis  is  to  use  the  damage  location  algorithm  to  determine  the 
location  of  the  damage.  Since  the  t)  row  norms  (defined  in  Section  5.3)  associated  with  this 

particular  FEM  are  of  the  same  order  of  magnitude,  either  the  subspace  rotation  direct 
method  or  the  angle  perturbation  method  can  be  used  to  locate  the  damage.  The  results  from 
the  application  of  the  subspace  rotation  direct  method  are  summarized  in  Figure  7.10.  The 
upper  left  plot  of  Figure  7.10  represents  the  exact  damage  location.  The  upper  right  plot 
corresponds  to  the  case  where  the  exact  (noise  free)  "damaged"  eigenvector  information  is 
provided  to  the  algorithm.  The  lower  left  and  right  plots  correspond  to  the  case  where  the 
exact  "damaged"  eigenvectors  have  been  corrupted  with  2.5%  and  5%  random  noise, 
respectively.  As  shown  in  Figure  7.10,  the  location  algorithm  is  able  to  exactly  locate  the 
damage  when  presented  with  noise  free  information.    It  is  also  seen  that  the  algorithm 
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experiences  a  degradation  in  detecting  damage  as  thie  error  in  the  eigenvectors  is  increased; 
however  the  damage  is  still  locatable. 
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Figure  7. 10      Fifty-Bay  Truss:  Damage  Location  Results  using  the  Subspace 
Rotation  Algorithm  with  the  Eigendata  of  the  First  Ten  Modes. 


7.3.3     Damage  Extent 

With  knowledge  of  the  location  of  damage  from  the  previous  analysis,  the  rank  of  the 
"true"  perturbation  to  the  stiffness  matrix,  AKd  can  be  found  by  adding  the  rank  of  the 
element  stiffness  matrix  of  the  damaged  members.  Hence,  the  rank  of  the  perturbation  to 
the  stiffness  matrix  due  to  damage  is  two  because  two  members  having  rank  one  element 
matrices  are  damaged.  From  Proposition  6. 1  and  6.2,  it  is  clear  that  only  experimental 
eigendata  from  two  modes  of  vibration  should  be  used  to  compute  the  extent  of  the  damage. 
In  the  noisy  cases,  the  two  modes  that  should  be  used  are  the  ones  with  dj's  that  most  cleanly 
demonstrate  the  damage  shown  in  Figure  7.10.    These  modes  can  be  determined  by 
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inspecting  the  individual  damage  vectors  di  associated  with  the  measured  modes.  In 
reviewing  the  individual  dj,  it  was  determined  that  modes  8  and  9  provide  the  best  insight 
into  the  state  of  the  damage.  The  results  of  applying  the  extent  algorithm  to  determine  the 
perturbations  to  the  stiffness  matrix  due  to  the  damage,  AK^,  are  shown  in  Figure  7.11.  The 
mesh  plots  are  3-dimensional  representations  of  the  perturbation  matrices.  The  rows  and 
columns  of  the  mesh  plots  correspond  to  the  rows  and  columns  of  the  perturbation  matrices. 
The  "height"  of  each  peak  represents  the  magnitude  of  the  perturbation  made  to  each  matrix 
element.  In  Figure  7.11,  the  upper  left  mesh  plot  represents  the  exact  damage.  The  upper 
right  plot  corresponds  to  the  case  where  exact  eigendata  are  used  in  computing  the  extent  of 
damage.  Note  that  with  only  two  noise  free  eigenvalues  /  eigenvectors  (modes  8  and  9) ,  the 
algorithm  is  able  to  reproduce  the  exact  damage.  As  stated  and  proved  in  Proposition  6.2, 
the  exact  damage  can  be  computed  using  any  two  noise  free  modes  that  have  a  corresponding 
B  (defined  in  Section  6.3.2)  of  rank  two.  The  left  and  right  mesh  plots  of  the  second  row 
of  Figure  7. 1 1  correspond  to  the  cases  where  the  eigenvectors  have  been  corrupted  with  2.5% 
and  5%  random  noise,  respectively.  The  algorithm  again  demonstrates  good  performance 
when  faced  with  noisy  eigendata  (middle  plots).  The  results  of  applying  the  eigenvector 
filtering  algorithm  described  in  Section  6.6.2  to  predict  the  extent  using  modes  8  and  9  are 
shown  in  the  left  and  right  lower  plots  of  Figure  7.11.  It  is  obvious  that  the  filtering  process 
greatly  enhances  the  accuracy  of  the  extent  algorithm.  Figure  7.12  displays  the  MRPT's 
assessment  of  the  damage  extent,  for  the  noisy  eigenvector  cases,  when  the  eigendata  of  the 
first  ten  modes  of  vibration  are  used.  A  comparison  of  the  results  in  Figure  7.12  and  those 
in  the  second  row  of  Figure  7.11  shows  that  the  above  procedure  which  resulted  in  the  use 
of  only  modes  8  and  9  to  assess  the  extent  provided  cleaner  results.  In  the  mesh  plots  of 
Figure  7.12,  additional  numerical  values  appear  at  the  "undamaged"  DOFs  which  can  be 
attributed  to  the  noise  added  by  the  addition  of  the  extra  modes.  The  resultant  perturbation 
matrices  when  using  all  ten  modes  are  of  rank  10  which  is  not  consistent  with  the  "actual" 
damage  perturbation  which  is  of  rank  2.  Furthermore,  this  "ten  modes  procedure"  is  less 
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Exact  Damage 


Modes  8  &  9 


Modes  8  &  9  +  2.5%  error 


Modes  8  &  9  +  5%  error 


Modes  8  &  9  +  2.5%  error 
Filtered 


Modes  8  &  9  +  5%  error 
Filtered 


Figure  7.11      Fifty -Bay  Truss:  Damage  Extent  Results  using  the  MRPT 
with  the  Eigendata  of  Modes  8  and  9. 


efficient  since  it  requires  the  inversion  of  a  1 0-by- 1 0  matrix,  contrasted  to  the  2-by-2  inverse 
used  in  the  "two  modes  procedure."  The  percentage  errors  of  the  predicted  damage  extent 
with  respect  to  the  exact  stiffness  damage  for  all  studied  cases,  with  and  without  filtering, 
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are  listed  in  Table  7.1.  It  can  be  seen  from  the  table,  the  results  obtained  when  the  eigenvector 
filtering  process  is  used  are  much  more  accurate.  The  results  from  the  "two  modes 
procedure"  are,  in  all  cases,  but  one,  better  than  the  "ten  modes  procedure". 


Modes  1  to  10  +  2.5%  error 


Modes  1  to  10  +  5%  error 


Figure  7.12      Fifty-Bay  Truss:  Damage  Extent  Results  using  the  MRPT 
With  the  Eigendata  of  the  First  Ten  Modes. 


Modes  8  &  9 


Modes  1  to  10 


Figure  7.13      Fifty-Bay  Truss:  Damage  Extent  Results  using  Baruch's  Algorithm. 


For  the  sake  of  comparison,  results  from  the  application  of  Baruch  and  Bar  Itzhack's 
algorithm  (1978)  are  provided  in  Figure  7.13.  These  results  were  obtained  using  (i)  noise 
free  eigendata  from  modes  8  and  9,  and  (ii)  noise  free  eigendata  from  the  first  ten  modes. 
In  both  cases,  Baruch's  algorithm  fails  to  ascertain  the  extent  of  the  damage.  In  the  first  case, 
the  approach  completely  "smeared"  the  changes  throughout  the  entire  stiffness  matrix.  In 
the  second  case,  the  changes  are  more  localized  around  the  damaged  areas. 
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7.4    Experimental  Study:  The  NASA  8-bay  Truss 

7.4.1     Problem  Description 

The  eight  bay  hybrid-scaled  truss  used  in  this  investigation  is  part  of  a  series  of  structures 
designed  for  research  in  dynamic  scale  model  ground  testing  of  large  structures  at  the  NASA 
Langley  Research  Center.  Among  other  studies,  a  complete  analytical  and  experimental 
analysis  of  this  truss  was  performed  to  generate  a  realistic  testbed  for  structural  damage 
location/extent  algorithms  (Kashangaki  et  al.,  1992;  Kashangaki,  1992).  The  truss 
configuration  used  in  this  analysis  was  cantilevered  and  instrumented  with  ninety-six 
accelerometers  to  measure  all  three  translational  degrees  of  freedom  at  each  of  the  thirty-two 
unconstrained  nodes.  A  schematic  of  this  truss  is  shown  in  Figure  7. 14.  Figure  7.15  shows  a 
close-up  of  the  truss  lacing  pattern  along  with  strut  definitions  that  will  be  used  in  the 
remainder  of  this  study. 


Bay  One 


Bay  Eight 


Figure  7. 14      The  NASA  Eight-Bay  Hybrid-Scaled  Truss:  Damage  Cases. 
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Figure  7.15      The  NASA  8-Bay  Truss:  Lacing  Pattern. 

A  ninety-six  degree  of  freedom  undamped  finite  element  model  (FEM)  of  the  original 
"healthy"  truss  was  generated  using  MSC/NASTRAN.  In  this  FEM,  each  truss  strut  was 
modeled  as  a  rod  element.  Concentrated  masses  were  added  at  each  node  to  account  for  the 
joint  and  instrumentation  mass  properties  (Kashangaki,  1992).  The  mass  and  strut  properties 
of  the  truss  are  summarized  in  Tables  7.2  and  7.3. 


Table  7.2    Mass  Properties  of  the  Eight  Bay  Truss. 


Weight 
(Pounds) 

Total  Number  in 
Eight  Bay  Truss 

Total  Weight 
(Pounds) 

Node  Ball 

0.0128 

36 

2.6935 

Longeron  Strut 

0.0396 

68 

2.4524 

Diagonal  Strut 

0.0598 

41 

0.4608 

Joint  Assembly 

0.0170 

218 

3.7001 

Triax  +  Block 

0.0243 

32 

0.7760 

Total  Weight  of  the  Truss 

10.0776 
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Table  7.3    Strut  Properties  of  the  Eight  Bay  Truss. 


Longeron 

Diagonal 

Strut  Length 

19.685  in 

27.839  in 

Strut  Stiffness  (EA/L)eff 

13,040  lb/in 

9,013  lb/in 

Experimental  modal  analysis  of  the  truss  was  performed  for  the  "no  damage"  and 
sixteen  damage  cases.  In  the  testing,  the  excitation  source  was  provided  by  two  shakers. 
Modal  parameters  were  identified  from  the  measured  frequency  response  functions  using  the 
Polyreference  complex  exponential  technique.  For  each  case,  five  modes  of  vibration  were 
identified.  Each  measured  mode  consists  of  a  natural  frequency  and  its  corresponding  mode 
shape  with  measurements  at  all  ninety-six  FEM  degrees  of  freedom.  For  the  "no  damage" 
case,  the  first,  second,  fourth  and  fifth  mode  are  bending  modes;  and  the  third  mode  is  the 
first  torsional  mode.  In  the  damage  cases,  the  same  five  modes  were  measured;  however 
their  order  varies,  since  in  some  cases,  the  damage  could  cause  mode  switching.  It  should  be 
noted  that  during  the  testing  process,  two  accelerometers  at  node  fifteen  failed.  These 
defective  accelerometers  would  affect  the  component  of  the  measured  modes  associated 
with  DOFs  44  and  45. 

Data  for  fifteen  of  the  sixteen  damage  cases  tested  were  received  from  NASA.  The 
fifteen  damage  cases  that  were  received  are  identified  in  the  schematic  of  the  truss  shown  in 
Figure  7.14.  For  each  damage  case,  the  type  of  element  and  the  FEM  degrees  of  freedom 
affected  by  damage  are  shown  in  Table  7.4. 

In  damage  cases  A  to  N,  structural  damage  consisted  of  the  full  removal  of  one  strut 
from  the  truss.  Case  O  damage  consists  of  the  full  removal  of  two  struts.  In  Case  P,  one  of  the 
struts  was  buckled  to  illustrate  a  partial  damage  scenario.  Note  that  the  forthcoming  analysis 
of  the  fifteen  damage  cases  was  performed  with  no  a  priori  knowledge  of  the  actual  damage 
locations. 
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Table  7.4    NASA  8-Bay  Truss  Damage  Case  Definitions. 


Damage  Case 

Element  Type 

Damaged  DOFs 

A 

Longeron 

86 

C 

Longeron 

74,86 

D 

Upper  Diagonal 

73,  74,  88,  89 

E 

X-Batten 

73,76 

F 

Face  Diagonal 

76,78,79,81 

G 

Longeron 

65,77 

H 

Longeron 

62,74 

I 

Left  Diagonal 

68,  69,  74,  75 

J 

Longeron 

59,71 

K 

Right  Diagonal 

41,42,59,60 

L 

Longeron 

44,56 

M 

Longeron 

29,41 

N 

Longeron 

2,14 

O* 

Longeron  &  Diagonal 

62,  68,  69,  74,  75 

p- 

Z-Batten 

75,81 

*      Two  struts  removed 
**     Buckled  strut 


7.4.2     Refinement  of  the  Original  FEM 

The  original  FEM  provided  by  NASA  is  a  first  generation  analytical  model  that  was  not 
modified  to  match  the  experimentally  measured  modal  parameters.  It  was  found  that  the 
original  FEM  does  not  accurately  predict  the  dynamic  behavior  of  the  actual  structure.  As 
an  illustration  to  the  inaccuracy  of  the  original  FEM,  Table  7.5  shows  a  comparison  between 
the  analytically  computed  and  experimentally  measured  natural  frequencies.  Another 
commonly  used  criterion  to  judge  the  accuracy  of  the  original  analytical  FEM  is  to  compare 
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its  frequency  response  functions  to  those  measured  experimentally.  One  such  comparison  is 
shown  in  Figure  7.16.  It  is  clear  from  this  typical  frequency  response  comparison  that  the 
original  FEM  is  lacking  in  accuracy. 

Table  7.5    Comparison  of  Analytical  and  Experimental  Frequencies 


Natural 

Frequencies  (Hz) 

MODE# 

Analytical 

Experimental 

1 

13.9245 

13.8757 

2 

14.4407 

14.4783 

3 

46.7445 

48.4122 

4 

66.0067 

64.0332 

5 

71.1420 

67.4631 

Hence,  the  first  step  in  this  study  was  to  refine  the  original  FEM.  For  the  refinement 
process,  it  was  assumed  that  the  original  mass  matrix  is  an  accurate  representation  of  the 
truss'  mass  properties.  The  inaccuracy  of  the  original  FEM  was  believed  to  be  solely  due  to 
modeling  errors  in  the  stiffness  properties.  The  minimum  rank  perturbation  theory  discussed 
in  Chapter  6  (Secfion  6.3.2)  was  used  to  correct  the  original  stiffness  matrix.  In  order  to  get  a 
symmetric  updated  stiffness  matrix,  the  measured  mode  shapes  (eigenvectors)  were  mass 
orthogonalized  using  the  Optimum  Weighted  Orthogonalization  technique  (Baruch  1978). 

A  mesh  plot  of  the  changes  made  to  the  stiffness  matrix  (AK)  from  the  refinement 
process  is  shown  in  Figure  7.17.  As  can  be  seen  in  the  figure,  there  are  three  areas  of  major 
change  made  to  the  original  stiffness  matrix.  One  area  of  change  corresponds  to  the 
cantilever  end  of  the  truss.  This  change  was  expected  since  a  perfect  cantilever  condition  was 
assumed  in  generating  the  original  FEM.  Note  that  it  is  common  knowledge  that  perfect 
cantilever  conditions  cannot  be  produce  in  practice.  Other  changes  can  be  seen  in  the  middle 
of  the  mesh  plot.  These  changes  occur  at  and  around  DOFs  44  and  45  which  are  the  location 
of  the  bad  sensors.  Hence,  these  changes  were  attributed  to  the  effect  of  the  bad  sensors.  Yet, 
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another  area  of  changes  can  be  seen  near  the  free  end.  These  areas  coincide  with  the  locations 
of  the  two  shakers.  The  reason  for  the  perturbation  at  these  shaker  locations  is  that  their 
localized  stiffness  effects  were  not  modelled  in  the  original  FEM.  The  improvement  made  to 
the  original  analytical  FEM  from  the  MRPT  refinement  process  can  also  be  seen  in  Figure 
7.16.  It  is  apparent  from  comparing  the  frequency  response  functions  corresponding  to 
experimental  measurement,  original  FEM  and  the  refined  FEM  that  the  MRPT  provided  a 
great  deal  of  improvement  to  the  original  FEM.  The  amplitude  mismatch  between  the 
experimental  data  and  refined  FEM  prediction  can  be  attributed  to  the  fact  that  the  actual 
effects  of  the  damping  were  not  considered  in  the  analytical  models.  Note  that  the  actual 
damping  of  the  truss  is  small  and  negligible  but  it  is  nonzero. 


106 

105 
104 

103 

o 

M 

10' 
100 
10-1 

10-2 


1 1 

1 

1                          1                          1 

-  Experimental  Measurement 

—  Original  Analytical  Model 

ij 

! 

Hi 

M  '                               ::  !'  ^- 

..  Modified  Analytical  Model 

!  1 

'  \ 

1    ^ 

■A' 

\ 

1 

20  40  60  80 

Frequency  (Hz) 


100 


120 


140 


Figure  7.16      NASA  8-bay  Truss:  Typical  Frequency  Response  Comparison. 
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Figure  7.17      NASA  8-bay  Truss:  Perturbation  to  the  Original  Stiffness 
Matrix  that  Resulted  from  the  Refinement  Process. 


7.4.3     Damage  Location 

For  each  damage  case,  the  modal  parameters  of  the  five  measured  modes  were  used  to 
compute  a  cumulative  damage  location  vector  (CDLV)  to  determine  the  location  of  the 
structural  damage.  Since  the  values  of  ||  z>^  \\,  as  defined  in  Section  5.3,  are  of  the  same  orders 

of  magnitude  for  all  DOFs  and  all  five  measured  modes,  the  cumulative  damage  location 
vector  (CDLV)  could  be  computed  using  either  Eq.  (5.4.1)  or  Eq.  (5.4.2).  The  angle 
perturbation  CDLV  (Eq.  (5.4.2))  associated  to  damage  case  F  is  shown  in  Figure  (7.18).  The 
upper  left  comers  of  Figures  (7.19)  to  (7.32)  display  plots  of  the  cumulative  damage  location 
vectors  as  calculated  by  Eq.  (5.4.2)  for  all  cases  where  the  damage  was  successfully  located. 
In  these  plots,  the  CDLV  were  unit  normalized  and  their  elements  were  plotted  versus  the 
FEM  DOFs.  The  performance  of  the  location  algorithm  for  each  damage  case  is  summarized 
in  the  second  column  of  Table  7.7.  A  summary  of  the  procedures  used  in  the  interpretation  of 
the  CDLV  of  each  damage  case  is  as  follows: 
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(1)  Initially,  a  comparative  study  of  the  CDLV  for  all  cases  was  performed  to  deduce  the 
effects  of  the  bad  accelerometers.  A  fairly  large  numerical  component  at  DOF  44  was 
detected  in  most  of  the  fifteen  CDLV.  This  component  was  believed  to  be  due  to  the 
"bad"  sensor  located  at  DOF  44.  Hence,  for  all  cases,  it  was  decided  to  ignore  the 
component  of  the  CDLV  at  DOF  44. 

(2)  For  damage  cases  A,  C,  D,  G,  H,  I,  J,  K,  M,  N  and  O,  the  location  of  the  damage  was 
determined  by  simply  considering  the  DOFs  associated  with  the  substantially  larger 
numerical  components  of  the  CDLV's  as  being  the  damaged  DOFs.  The  smaller 
numerical  components  of  the  CDLV's  at  the  other  DOFs  can  be  attributed  to 
measurement  errors.  Notice  that  for  cases  D,  I,  J,  M  and  N,  the  component  of  the  CDLV 
at  DOF  44  is  fairly  small  and  hence  negligible.  It  should  be  also  noted  that  the  damage  in 
Case  A  affects  only  one  degree  of  freedom  (86)  since  the  strut  involved  is  a  longeron 
connected  to  the  cantilevered  end. 

(3)  For  damage  case  E,  the  DOFs  associated  with  the  CDLV  components  of  greater  order  of 
magnitude  are  42, 73, 75, 76  and  8 1 .  By  utilizing  the  connectivity  of  the  original  FEM, 
the  combinations  of  these  DOFs  that  are  physically  meaningful,  i.e.  bound  a  strut,  are 
(73, 76)  and  (75,  81).  Based  on  these  results,  it  was  decided  that  two  struts  connecting 
DOF  (73, 76)  and  (75, 81)  were  damaged.  However,  in  actuality  the  only  damaged  strut 
is  the  one  connecting  DOFs  (73,  76). 

(4)  For  damage  case  L,  by  ignoring  DOF  44  ( the  "bad"  sensor),  the  only  DOF  that  is  clearly 
damaged  is  56.  No  strut  with  DOF  56  is  connected  to  the  wall.  The  damage  of  a  strut  not 
connected  to  the  cantilevered  end  must  affect  at  least  two  DOFs .  Since  DOF  56  is  in  the 
y-direction,  it  was  deduced  that  the  most  probable  damaged  strut  is  a  longeron 
connected  at  one  end  to  DOF  56.  The  two  candidate  struts  are  the  ones  bounded  by 
DOFs  (44,  56)  and  (56,  68).   However,  the  component  of  the  CDLV  associated  with 
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DOF  68  is  small  and  is  in  no  way  affected  by  damage.  This  prompted  the  deduction  that 
the  strut  connecting  DOFs  44  and  56  is  the  damaged  member  in  Case  L. 

(5)  For  damage  case  P,  the  substantially  larger  components  of  the  CDLV  occur  at  DOFs  42, 
75  and  81  (DOF  44  ignored).  Again,  with  the  use  of  the  connectivity  of  the  original 
FEM,  it  is  deduced  that  the  only  combination  of  these  three  DOFs  that  are  connected  by  a 
strut  is  (75,  81).  These  DOFs  are  exactly  the  damaged  DOFs  of  Case  P. 

(6)  Damage  Case  F  is  the  only  case  where  the  location  algorithm  failed  to  locate  the  damage. 
The  cumulative  damage  vector  for  this  case  is  shown  in  Figure  7.18.  The  reason  the 
damage  was  not  located  in  this  case  is  that  it  involves  a  face  diagonal  strut.  Face 
diagonal  struts,  in  general,  affect  only  the  axial  modes  and  have  slight  or  no  effects  on 
the  bending  and  torsional  modes.  As  reported  earlier,  the  only  available  (measured) 
modes  are  bending  and  torsional  modes.  It  was  determined  analytically  that  the  first 
axial  mode  of  the  truss  occurs  at  the  sixth  mode. 

Kashangaki  et  al.  (1992)  performed  a  pure  analytical  study  to  compile  a  hst  of  the 
elements  of  the  8-bay  truss  that  have  a  substantial  contribution  to  the  total  strain  energy  for 
each  of  the  five  measured  modes.  The  present  study  confirms  that  the  damage  in  struts  that 
have  a  substantial  contribution  to  the  total  modal  strain  energy  for  the  measured  modes  is 
detectable.  However,  our  study  also  shows  that  damage  of  non-highly  strained  struts  for  the 
measured  modes  can  also  be  located  as  illustrated  in  cases  G  and  N.  Note  that  the  damage  in 
cases  E  and  F  involves  struts  with  low  strain  energy  for  the  measured  modes. 

7.4.4     Damage  Extent 

Because  the  truss  under  investigation  is  light  weight  and  very  flexible,  the  effect  of 
structural  damage  is  considered  (i)  negligible  on  the  overall  inertial  properties,  and  (ii)  a 
substantial  stiffness  loss.  Using  the  extent  algorithm  discussed  in  Section  6.3.2,  the 
perturbation  to  the  stiffness  matrix  (AKj)  due  to  the  damage  is  calculated  for  each  damage 
case  where  the  damage  was  successfully  located.    To  insure  the  symmetry  of  AK^^,  the 
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eigenvectors  used  in  these  calculations  were  first  mass  orthogonalized  using  the  Optimum 
Weighted  Orthogonalization  technique  (Baruch  and  Bar  Itzhack,  1978). 

7.4.4.1  The  Brute  Force  Method 

The  upper  right  comers  of  Figures  7.19  to  7.32  display  mesh  plots  of  the  calculated 
AKd's  using  the  modal  parameters  of  all  five  measured  modes.  From  these  mesh  plots,  it  is 
clear  that  on  the  average  the  extent  algorithm  failed  to  concentrated  the  major  changes  at  the 
damaged  DOFs.  This  unsatisfactory  performance  of  the  extent  algorithm  can  be  attributed  to 
the  measurement  error  of  the  defective  accelerometers.  In  most  cases,  the  major  changes  are 
concentrated  at  DOF  44,  which  corresponds  to  the  location  of  a  defective  accelerometer.  In 
addition,  for  any  given  damage  case,  as  discussed  in  Section  6.6.1,  only  a  subset  of  the 
measured  modes  are  highly  affected  by  damage.  Using  modes  that  are  not  highly  affected  by 
damage  only  introduces  their  associated  measurement  noise  in  the  extent  calculation. 

7.4.4.2  The  Damage  Consistent  Method 

As  discussed  in  Section  6.6. 1  and  illustrated  in  the  previous  example,  for  each  damage 
case,  the  performance  of  the  extent  algorithm  can  be  improved  by  only  using  as  many  modes 
as  the  rank  of  the  "actual"  AK^.  The  rank  of  the  "actual"  AIQ  can  be  estimated  by  adding  the 
rank  of  all  element  stiffness  matrices  which  connect  damaged  DOFs.  Since  the  elemental 
stiffness  matrix  of  the  truss  under  investigation  is  rank  one,  the  rank  of  the  "actual"  AKd,  for 
each  damage  case,  is  equal  to  the  number  of  damaged  struts.  Hence,  only  modal  parameters 
from  one  mode  are  needed  to  compute  the  damage  extent  of  all  cases  featuring  one  damaged 
strut.  For  Case  O,  since  two  struts  are  damaged,  data  from  two  modes  are  needed.  For  each 
case,  the  modes  that  should  be  used  are  the  ones  that  most  cleanly  demonstrate  the  state  of 
damage  of  that  case.  These  modes  can  be  simply  determined  by  inspecting  the  individual 
damage  location  vector  dj.  The  modes  that  provide  the  best  insight  into  the  state  of  the 
damage  for  each  damage  case  are  reported  in  the  third  column  of  Table  7.7.  The  unit 
normalized  damage  location  vectors  associated  with  the  "best"  modes  for  all  cases  are 
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plotted  versus  the  FEM  DOFs  in  the  middle  left  of  Figures  7. 19  to  7.32.  The  mesh  plots  of  the 
AKd's  computed  using  only  the  "best"  mode(s)  for  all  damage  cases  are  shown  in  the  middle 
right  comers  of  Figures  7. 19  to  7.32.  It  is  clear  that  the  performance  of  the  extent  algorithm 
for  all  cases  has  been  greatly  improved  in  comparison  to  its  performance  when  all  five 
measured  modes  were  used.  The  estimated  effective  stiffness  loss  (lb/in)  for  all  damage 
cases,  computed  from  this  process,  are  reported  in  the  fourth  column  of  Table  7.7.  For  each 
damage  case,  the  effective  stiffness  loss  was  calculated  by  averaging  the  stiffness  changes 
(AKd )  associated  with  the  damaged  DOFs.  In  the  original  FEM,  the  longerons,  X-battens 
and  Z-battens  have  an  effective  stiffness  of  13,040  lb/in,  while  the  effective  stiffness  of  the 
diagonal  struts  is  9,013  lb/in  (Kashangaki,  1992).  The  removal  of  any  strut  results  in  the 
complete  loss  of  the  stiffness  of  that  strut.  The  fifth  column  of  Table  7.7  reports  the 
percentage  error  of  the  estimated  stiffness  losses  with  respect  to  the  original  FEM.  The 
percentage  error  for  Case  P  (buckled  strut)  was  not  computed  since  the  damage  extent  is 
unknown. 

7.4.4.3    Application  of  the  Eigenvector  Filtering  Algorithm 

The  damage  extent  assessment  can  be  further  improved  by  using  the  eigenvector 
filtering  algorithm  discussed  in  Section  5.4.2.  For  consistency,  the  mode(s)  used  in  the 
filtering  process  are  the  "best"  mode(s)  as  determined  above  and  reported  in  the  third  column 
of  Table  7.7.  For  each  case,  the  "undamaged"  components  of  the  damage  vectors  d,  are  first 
set  to  zero.  The  components  of  di  associated  to  the  "damaged"  DOFs  are  then  constrained  to 
be  consistent  with  the  actual  effect  of  the  damage  as  dictated  by  the  element  stiffness  matrix 
of  the  damage  struts.  The  damage  in  X-battens,  longerons  or  Z-battens  affects  two  DOFs 
(one  DOF  at  each  of  the  two  nodes  of  the  strut).  In  order  for  the  damage  to  be  finite  element 
consistent,  the  component  associated  with  these  two  DOFs  should  be  equal  in  magnitude  and 
of  opposite  sign.  Damage  in  upper,  face,  left  or  right  diagonal  struts  affects  four  DOFs  (two 
at  each  node  of  the  strut) .  These  four  DOFs  can  be  classified  into  two  sets  of  two  DOFs  in  the 
same  direction  (x,  y  or  z).  The  DOFs  of  a  given  set  are,  in  general,  equal  in  magnitude  and  of 
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opposite  signs.  In  the  8-bay  truss,  all  upper,  face,  left  or  right  diagonals  are  either  at  a  45°  or 
-45°  angle  with  respect  to  their  in-plane  global  coordinate  axis.  This  constrains  the  DOFs  of 
the  two  sets  to  have  the  same  magnitude  in  order  to  be  consistent  with  the  raM.  The  DOFs 
associated  with  a  given  node  have  the  same  sign  when  the  angle  is  45°;  they  are  of  opposite 
sign  when  the  angle  is  ^5°.  For  the  cases  where  a  single  strut  is  damaged  the  magnitude  of  a 
DOF  could  be  set  arbitrarily.  The  reason  is  that  eigenvectors  are  unique  in  a  relative  sense; 
i.e.,  if  Ui  is  an  eigenvector  associated  to  a  given  mode  then  ajii  is  also  an  eigenvector 
associated  to  that  same  mode,  a  being  any  scalar.  Thus,  the  components  of  the  damage  vector 
are  also  a  function  of  the  scalar  a.  These  components  can  be  set  to  any  values  by  varying  a. 
For  the  single  member  damage  cases  under  study,  this  process  is  sunmiarized  in  Table  7.6. 


Table  7.6    Summary  of  the  Filtering  Process  for  Single  Member  Damage  Cases. 


Node  1 

Node  2 

DOF 

X 

Y 

z 

X 

Y 

z 

X-Batten 

P 

0 

0 

-P 

0 

0 

Longeron 

0 

P 

0 

0 

-p 

0 

Z-Batten 

0 

0 

p 

0 

0 

-p 

Case  D  M5°) 

p 

-p 

0 

-p 

p 

0 

Case  I  (45°) 

p 

0 

p 

-p 

0 

-p 

Case  K  M5°) 

p 

0 

-p 

-p 

0 

p 

where  P  is  any  scalar 


Note  that  cases  A,  C,  E,  G,  H,  J,  L,  M,  N  and  P  involve  the  damage  of  either  an  X-batten, 
longeron  or  Z-batten.  As  discussed  earlier,  all  DOFs  of  the  filtered  damage  vector  d^ , 
associated  to  a  given  damage  member  have  the  same  magnitude,  p.  Experience  gained  by 
using  the  eigenvector  filtering  algorithm  indicates  that  in  the  multiple  member  damage 
scenario  (case  O),  the  ratio  between  the  P's  from  different  damaged  struts  contains  important 
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information  about  the  relative  damage  extent  among  the  damaged  struts.  Hence,  in  case  O,  a 
given  P  associated  with  a  given  damaged  strut  cannot  be  chosen  arbitrary.  This  P  should  be 
estimated  by  averaging  the  components  of  the  unfiltered  damage  vector  d^  corresponding  to 

the  damaged  strut  in  question.  The  absolute  value  of  the  unit  normalized  filtered  damage 
location  vectors  associated  with  all  cases  are  plotted  versus  the  FEM  DOFs  in  the  bottom  left 
corners  of  Figures  7. 19  to  7.32.  Now  that  the  "filtered"  damage  vectors  are  generated,  the 
next  step  is  to  compute  the  corresponding  filtered  damage  eigenvector  using  Eq.  5.4.3.  The 
process  of  computing  the  perturbation  due  to  the  damage  using  the  filtered  eigenvectors  is  as 
in  the  previous  section.  For  all  damage  cases,  the  mesh  plots  of  the  AK^'s  computed  using 
the  "filtered"  mode(s)  are  displayed  in  the  bottom  right  comers  of  Figures  7. 19  to  7.32.  The 
effective  stiffness  loss  computed  using  the  filtering  process,  for  each  damage  case,  is  » 
reported  in  the  sixth  column  of  Table  7.7.  The  seventh  column  of  Table  7.7  reports  the  ' 
percentage  error  of  the  estimated  stiffness  losses  with  respect  to  the  original  raM.  A 
comparison  of  these  results  to  the  one  acquired  from  the  previous  section  shows  that  the 
eigenvector  filtering  algorithm  greatly  enhances  the  extent  estimate. 


Damage  Location:  Modes  1-5 


Figure  7.18      NASA  8-bay  Truss:  Cumulative  Damage  Vector 
Associated  to  Case  F. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 


Damage  Location  -  Mode  1 


Damage  Extent  -  Mode  1 


Damage  Location  -  Mode  1  Filtered 


Damage  Extent  -  Mode  1  Filtered 


Figure  7. 19      NASA  8-bay  Truss:  Damage  Assessment  of  Case  A. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 


Damage  Location  -  Mode  1 


Damage  Extent  -  Mode  1 


Damage  Location  -  Mode  1  Filtered  Damage  Extent  -  Mode  1  Filtered 


Figure  7.20      NASA  8-bay  Truss:  Damage  Assessment  of  Case  C. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 


Damage  Location  -  Mode  3 


Damage  Extent  -  Mode  3 


Damage  Location  -  Mode  3  Filtered 
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Damage  Extent  -  Mode  3  Filtered 


Figure  7.21      NASA  8-bay  Truss:  Damage  Assessment  of  Case  D. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 


Damage  Location  -  Mode  4 


Damage  Extent  -  Mode  4 


Damage  Location  -  Mode  4  Filtered 
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Figure  7.22      NASA  8-bay  Truss:  Damage  Assessment  of  Case  E. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 


Damage  Location  -  Mode  1 


Damage  Extent  -  Mode  1 


Damage  Location  -  Mode  1  Filtered 


Damage  Extent  -  Mode  1  Filtered 


Figure  7.23      NASA  8-bay  Truss:  Damage  Assessment  of  Case  G. 
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Damage  Location  -  Modes  1-5  Damage  Extent  -  Modes  1-5 


8)0., 


44 

iHlMiflibilk  iUi 

62 

\Mm 

11 

74        - 

0 


50 
DOF 


Damage  Location  -  Mode  4 
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Damage  Location  -  Mode  4  Filtered  Damage  Extent  -  Mode  4  Filtered 


Figure  7.24      NASA  8-bay  Truss:  Damage  Assessment  of  Case  H. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 


Damage  Location  -  Mode  4 


Damage  Extent  -  Mode  4 


Damage  Location  -  Mode  4  Filtered 


Damage  Extent  -  Mode  4  Filtered 


Figure  7.25      NASA  8-bay  Truss:  Damage  Assessment  of  Case  L 
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Figure  7.26      NASA  8-bay  Truss:  Damage  Assessment  of  Case  J. 
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Damage  Location  -  Modes  1-5 
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Damage  Extent  -  Mode  3 


Damage  Location  -  Mode  3  Filtered 


Damage  Extent  -  Mode  3  Filtered 


Figure  7.27      NASA  8-bay  Truss:  Damage  Assessment  of  Case  K. 
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Damage  Location  -  Modes  1-5 
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Figure  7.28      NASA  8-bay  Truss:  Damage  Assessment  of  Case  L. 
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Damage  Location  -  Modes  1-5 
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Damage  Extent  -  Modes  1-5 
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Damage  Location  -  Mode  4  Filtered  Damage  Extent  -  Mode  4  Filtered 


Figure  7.29      NASA  8-bay  Truss:  Damage  Assessment  of  Case  M. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 
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Figure  7.30      NASA  8-bay  Truss:  Damage  Assessment  of  Case  N. 
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Figure  7.3 1      NASA  8-bay  Truss:  Damage  Assessment  of  Case  O. 
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Damage  Location  -  Modes  1-5 


Damage  Extent  -  Modes  1-5 
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Figure  7.32      NASA  8-bay  Truss:  Damage  Assessment  of  Case  P. 
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Table  7.7    NASA  8-Bay  Truss:  Summary  of  the  Damage  Assessment  Results. 


Without  Filtering 

With  Filtering 

Damage 
Case 

Location 
Performance 

Best 
mode 

Extent 
(lb/in) 

%  error 

Extent 
(lb/in) 

%  error 

A 

++ 

1 

11,896 

8.77 

12,865 

1.34 

C 

++ 

1 

11,901 

8.73 

12,866 

1.34 

D 

++ 

3 

9,789 

8.61 

8,806 

2.30 

E 

o 

4 

8,094 

37.90 

13,454 

3.17 

F 

- 

Location  not  detected 

G 

++ 

1 

10,859 

16.70 

12,971 

0.53 

H 

++ 

4 

12,601 

3.37 

13,154 

0.87 

I 

++ 

4 

7,568 

16.00 

10,220 

13.39 

J 

++ 

4 

12,015 

7.86 

13,176 

1.04 

K 

++ 

3 

6,955 

22.80 

8,851 

1.79 

L 

+ 

4 

8,421 

35.40 

12,859 

1.39 

M 

++ 

4 

12,009 

7.91 

12,882 

1.21 

N 

++ 

5 

10,037 

23.00 

13,169 

0.99 

O 

++ 

3,5 

(L)  12,347 
(D)  8,451 

5.31 
6.78 

(L)  12,807 
(D)  8,878 

1.79 
1.50 

P 

+ 

4 

8,815 

not  def. 

13,461 

not  def. 

++   Damage  clearly  located 

+      Damage  located  with  further  analysis 

-      Damage  not  located 

o      location  narrowed  to  within  two  members 

(L)  Longeron  (D)  diagonal 
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7.5    Experimental  Study:  Mass  Loaded  Cantilevered  Beam 

7.5.1     Problem  Description 

This  experiment  was  designed  to  illustrate  the  scenario  of  damage  in  mass  properties. 
The  structure  used  in  this  investigation  is  a  cantilevered  beam  loaded  with  a  non-structural 
mass.  A  schematic  of  this  beam  is  shown  in  Figure  7.33.  The  dimensions  and  properties  of 
the  structure  are  summarized  in  Table  7.8. 


Discrete  Mass 


1 


1  2 


Figure  7.33      The  Mass  Loaded  Cantilevered  Beam. 


Table  7.8    Mass  Loaded  Cantilevered  Beam  Properties. 


Beam  Length  -  .86  m 

Beam  Mass/Length  -  1 .246  kg/m 

Beam  Moment  of  Inertia  -  1.458x10""^  m'^ 

Beam  Young's  Modulus  -  69  GPa 

Discrete  Mass  -  0.7938  kg 

Discrete  Mass  Moment  of  Inertia  -  1.1x10"^  kg-m^ 


An  undamped  FEM  of  the  mass  loaded  beam  was  constructed  using  beam  elements  in 
conjunction  with  the  properties  of  the  non-structural  mass.  The  beam  element  has  two 
degrees  of  freedom  (DOF)  at  each  node:  bending  and  rotation.  A  sixteen  DOF  undamped 
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FEM  was  generated  using  the  eight  equal  length  element  discretization  shown  in  Figure  7.33. 
The  effect  of  the  non  structural  mass  was  considered  non-stiffening  and  concentrated  at  node 
3.  This  mass  was  modelled  by  adding  its  mass  and  moment  of  inertia  to  node  3  bending  and 
rotation  DOFs,  respectively.  The  mass  loaded  beam,  as  described  above,  was  considered  the 
"healthy"  configuration.  Structural  damage  consisted  of  the  removal  of  the  non-structural 
mass  from  the  beam. 

Experimental  modal  analysis  of  the  beam  was  performed  on  the  "healthy  and 
"damaged"  beam  configurations.  Modal  parameters  were  identified  using  frequency 
domain  techniques  and  the  Rational  Fraction  Least  Square  single  degree  of  freedom  curve 
fitting  algorithm.  The  excitation  source  used  was  an  impact  hammer  and  the  driving  point 
measurement  was  an  accelerometer  mounted  at  the  free  end  of  the  beam.  For  each 
configuration  (healthy  and  damaged),  four  modes  of  vibration  were  measured.  Each  mode 
consisted  of  a  natural  frequency  and  its  corresponding  mode  shape  with  measurements  at 
only  the  eight  FEM  bending  degrees  of  freedom. 

7.5.2  Analytical  and  Experimental  Models  Dimension  Correlation 

The  number  of  measured  eigenvector  components  (8)  is  less  than  the  number  of  DOFs  in 
the  FEM  (16).  In  fact,  only  the  bending  DOFs  of  the  beam  were  measured  experimentally. 
As  discussed  in  Chapter  2,  two  approaches  are  available  to  correlate  these  dimensions:  (i) 
expansion  of  the  measured  eigenvectors  or  (ii)  reduction  of  the  FEM.  It  was  found  that  a 
FEM  reduction  is  better  suited  for  this  application.  Thus,  the  FEM  was  reduced  using  the 
improved  reduction  system  (IRS)  method  (O'Callahan,  1989). 

7.5.3  Refinement  of  the  Original  FEM 

Table  7.9  shows  a  comparison  between  the  analytically  computed  and  experimentally 
measured  natural  frequencies.  From  this  comparison,  it  is  clear  that  the  original  reduced 
FEM  does  not  accurately  predict  the  dynamic  behavior  of  the  "healthy"  beam.  The  first  step 
in  this  study  was  to  refine  the  original  FEM.  For  the  refinement  process,  it  was  assumed  that 
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the  original  mass  matrix  is  an  accurate  representation  of  the  structure's  mass  properties.  The 
inaccuracy  of  the  original  FEM  was  believed  to  be  solely  due  to  modeling  errors  in  the 
stiffness  properties.  The  algorithm  discussed  in  Section  6.3.2  was  used  to  correct  the  original 
reduced  stiffness  matrix.  In  order  to  get  a  symmetric  updated  stiffness  matrix,  the  measured 
mode  shapes  (eigenvectors)  were  mass  orthogonalized  using  the  Optimum  Weighted 
Orthogonalization  technique  (Baruch  and  Bar  Itzhack,  1978). 


Table  7.9  Analytical  and  Experimental  Frequencies 

of  the  "Healthy"  Structure. 


Mode# 

Analytical 
Frequency  (Hz) 

Experimental 
Frequency  (Hz) 

1 

4.5 

4.6 

2 

40.3 

41.6 

3 

96.8 

99.8 

4 

195.9 

205.9 

7.5.4       Damage  Location 

The  next  step  of  this  analysis  is  to  determine  the  location  of  the  structural  damage.  Using 
the  refined  FEM,  and  the  modal  parameters  of  the  four  modes  of  vibration  measured  from  the 
"damaged"  beam,  a  cumulative  damage  location  vector  was  calculated  (Chapter  5).  Since 
the  values  of  ||  ^  ||,  as  defined  in  Section  5.3,  are  of  different  orders  of  magnitude  the 

cumulative  damage  location  vector  should  be  computed  using  Eq.  (5.4.2).  The  upper  left 
comer  of  Figure  7.34  displays  the  plot  of  the  unit  normalized  cumulative  damage  location 
vector  as  calculated  by  Eq.  (5.4.2).  From  this  plot,  it  is  clear  that  DOF  3  has  been  affected  by 
damage.  This  is  exactly  the  bending  DOF  where  the  non-structural  mass  was  mounted.  The 
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small  numerical  elements  at  all  other  DOFs  can  be  attributed  to  experimental  measurement 
noise. 


5P      0.5 

e 


Damage  Location  -  Modes  1-4 


a   0.5 


Damage  Location  -  Mode  4 


4  6 

DOF 


4  6 

DOF 


Damage  Extent  -  Mode  4 


Figure  7.34      Mass  Loaded  Cantilevered  Beam:  Damage  Assessment 


7.5.5     Damage  Extent 

The  final  step  of  this  analysis  is  to  determine  the  extent  of  structural  damage.  With 
knowledge  of  the  damage  location,  the  rank  of  the  "true"  mass  perturbation  matrix,  AMj,  is 
one  since  there  is  only  one  DOF  affected  by  the  damage  and  that  DOF  is  not  connected  to  the 
cantilevered  end.  Hence,  in  order  to  compute  a  rank  one  AM^  only  one  mode  of  vibration  is 
needed.  As  discussed  earlier,  the  mode  that  should  be  used  is  the  one  that  most  cleanly 
demonstrates  the  damage  detected  by  the  cumulative  damage  location  vector.  The  damage 
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vector  d4  was  determined  to  provide  the  best  insight  into  the  state  of  damage.  The  damage 
vector,  d4  is  shown  in  the  upper  right  comer  of  Figure  7.34.  The  calculated  AMj  using  mode 
4  data  is  shown  in  the  lower  plot  of  Figure  7.34.  It  is  clear  that  the  extent  calculation  has 
concentrated  the  major  changes  at  the  DOF  affected  by  the  damage.  From  the  extent 
calculation,  the  mass  loss  was  estimated  to  be  0.7947  kilograms,  which  is  within  1 . 1  %  of  the 
actual  mass  loss. 

7.6    Fifty-Bay  Two-Dimensional  Truss:  Proportionally  Damped  FEM 

7.6.1     Problem  Description 

The  structure  used  in  this  investigation  is  the  same  as  the  one  used  in  Section  7.3.  For 
convenience,  it  is  shown  again  in  Figure  7.35.  This  example  is  used  to  illustrate  the 
characteristics  of  the  MRPT  when  dealing  with  proportionally  damped  systems  in  which  two 
of  the  three  property  matrices  are  affected  by  structural  damage.  The  analytical  model  of  the 
truss  is  a  201  DOF  FEM  which  was  generated  using  rod  elements.  The  damping  of  the 
"healthy"  model  is  proportional  and  equal  to  lxlO~'  times  the  "healthy"  mass  matrix  plus 
5x10"^  times  the  "healthy"  stiffness  matrix.  As  in  the  problem  of  Section  7.3,  damage  is 
simulated  by  reducing  the  Young's  modulus  of  two  members:  the  upper  longeron  of  bay  three 
and  the  lower  longeron  of  bay  forty.  The  modulus  of  elasticity  of  the  upper  longeron  of  the 
third  bay  is  reduced  from  E=29xl0^  psi  to  1x10^  psi.  The  bay  forty's  lower  longeron  is 
subjected  to  a  complete  loss  of  stiffness  (Young's  modulus  equal  to  zero).  The  damping  of 
the  "damaged"  model  is  also  proportional  and  equal  to  lxl0~^  times  the  "healthy"  mass 
matrix  plus  5x1 0~^  times  the  "damaged"  stiffness  matrix.  Thus,  the  simulated  damage  is 
only  affecting  the  stiffness  and  damping  properties  of  the  structure.  Note  that  the  damage  in 
the  damping  properties  is  proportional  to  the  damage  in  the  stiffness  properties.  For  our 
damage  assessment  analysis,  it  is  assumed  that  only  the  first  ten  "damaged"  modes  are 
available.  This  problem  is  investigated  for  three  different  scenarios.  Each  scenario 
corresponds  to  a  different  level  of  random  noise  added  to  the  "damaged"  eigenvectors.  In 
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practice  the  noise  in  eigenvector  information  could  be  due  to  both  measurement  and/or 
expansion  errors. 


Bay  1 


// 


// 


Bay  50 


E  =  2.0x101'  (29xl06psi)  p  =  7.833xl03  kg/m^  (0.283  Ib/in^) 

A  =  6.452X  10-4  m2  (1.0  in2)  L  =  0.1270  m  (5.0  in) 

Figure  7.35      50-Bay  2-Dimensional  Truss. 


7.6.2     Damage  Location 

Cumulative  damage  location  vectors,  as  defined  in  Chapter  5,  are  first  computed  for  all 
three  scenarios  using  the  modal  properties  of  the  ten  "damaged"  modes.  In  this  problem,  all 
rows  of  matrix  Z ^  (defined  in  Chapter  5)  are  of  the  same  order  of  magnitude;  hence  either 

Eq.  (5.4. 1)  or  Eq.  (5.4.2)  can  be  used  to  compute  the  cumulative  damage  vectors.  The  upper 
left  plot  of  Figure  7.36  represents  the  exact  damage  computed  from  the  exact  damage 
perturbation  matrices  (AIQ  and  ADj).  The  upper  right  plot  corresponds  to  the  case  where  the 
exact  "damaged"  eigenvector  information  is  provided  to  the  subspace  rotation  damage 
location  algorithm.  The  lower  left  and  right  plots  correspond  to  the  cases  where  the  exact 
"damaged"  eigenvectors  have  been  corrupted  with  2.5%  and  5%  random  noise,  respectively. 
As  shown  in  Figure  7.36,  the  location  algorithm  is  able  to  exactly  locate  the  damage  when 
presented  with  noise  free  information.  Although  not  as  clean,  the  damage  can  still  be  clearly 
located  in  the  noisy  eigenvectors  cases. 
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Figure  7.36      50-Bay  2-Dimensional  Truss:  Damage  Location. 


7.6.3     Damage  Extent 


The  damage  location  assessed  in  the  previous  section  can  be  used  in  conjunction  with  the 
truss  finite  element  connectivity  to  determine  the  damaged  truss  struts.  The  rank  of  the 
"true"  perturbation  matrix,  AKd  (or  ADj) ,  can  be  found  by  adding  the  rank  of  the  element 
stiffness  (or  damping)  matrix  of  the  damaged  struts.  Hence,  the  rank  of  the  perturbation  to 
the  stiffness  matrix  due  to  damage  is  two  because  two  struts  of  rank  one  element  matrices  are 
damaged.  Note  that  the  element  stiffness  matrix  of  a  strut  is  one  since  it  is  modelled  as  a  rod 
element.  Because  the  damping  of  the  structure  is  proportional  and  the  mass  properties  are 
unaffected  by  the  damage,  it  is  deduced  that  the  rank  of  the  perturbation  to  the  damping 
matrix  is  the  same  as  the  rank  of  the  perturbation  to  the  stiffness  matrix.  From  propositions 
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6. 1  and  6.2,  it  is  clear  that  only  experimental  data  from  two  modes  of  vibration  are  needed  to 
compute  the  extent  of  the  damage.  In  the  noisy  situations,  the  two  modes  that  should  be  used 
are  the  ones  that  most  cleanly  demonstrate  the  damage  shown  in  Figure  7.36.  These  modes 
can  be  determined  by  inspecting  the  individual  damage  vectors.  An  inspection  of  the 
individual  damage  vectors  associated  with  each  "noisy"  eigenvector  suggests  that  modes  8 
and  9  provide  the  best  insight  into  the  state  of  the  damage.  The  results  of  applying  the  MRPT 
(  Section  6.4.1,  Eq.  (6.4.6))  to  determine  the  perturbations  to  the  stiffness  matrix  due  to  the 
damage  are  shown  in  Figure  7.37.  Clearly,  from  judging  the  upper  right  mesh  plot  of  Figure 
7.37,  the  MRPT  is  able  to  reproduce  the  exact  damage  with  only  two  noise  free  modes.  The 
algorithm  demonstrates  good  performance  when  faced  with  noisy  eigendata  (lower  plots). 
The  percentage  errors  with  respect  to  the  exact  stiffness  damage  for  all  studied  cases  are 
listed  in  Table  7.10.  The  perturbations  to  the  damping  matrix,  ADj,  due  to  the  damage  are 
estimated  by  the  extent  algorithm  (Section  6.4. 1 ,  Eq.  (6.4.4))  with  exactly  the  same  accuracy 
as  for  AKd.  In  fact,  the  unsealed  mesh  plots  of  the  computed  ADd's  are  the  same  as  the  ones 
shown  in  Figure  7.37.  This  would  be  expected  since,  as  reported  earlier,  the  damage  in  the 
damping  properties  is  proportional  to  the  damage  in  the  stiffness  properties. 
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Table  7. 10     50-Bay  2-Dimensional  Truss:  Summary  of  Percentage 
Error  With  Respect  to  the  Exact  Damage. 


Percentage  Error  with  respect  to 
exact  stiffness  (or  damping) 

Eigenvectors 
Error 

Upper  Longeron  Of 
Bay  Three 

Lower  Longeron  of 
Bay  Forty 

0.0% 

0.00 

0.00 

2.5% 

6.60 

7.42 

5.0% 

22.60 

15.30 
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Exact  Damage 


Modes  8  &  9 


Modes  8  &  9  +  2.5%  error 


Modes  8  &  9  +  5%  error 


Figure  7.37      50-Bay  2-Dimensional  Truss:  Damage  Extent. 


7.7    Eight-Bay  Two-Dimensional  Mass-Loaded  Cantilevered  Truss 


7.7.1     Problem  Description 

The  structure  under  investigation  is  the  eight-bay  two-dimensional  mass-loaded 
cantilevered  truss  shown  in  Figure  7.38.  This  structure  was  designed  (Ricles,  1991)  to 
emulate  typical  properties  of  space  structures:  low  frequency  modes  with  large 
non-structural  mass.  The  geometric  and  material  properties  of  the  truss  are  given  in  Figure 
7.38.  The  truss  consists  of  40  struts  and  16  nodes.  Fourteen  of  the  nodes  are  loaded  with 
concentrated  non- structural  mass  of  magnitude  2.3  Ibs-sec^/in  and  the  remaining  two  nodes 
have  large  lumped  masses  of  magnitude  47.4  Ibs-sec^/in.  Each  truss  strut  was  modeled  as  a 
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rod  element  with  negligible  mass.  The  finite  element  model  has  32  translational  DOFs  ( 2 
DOFs  per  node).  The  truss,  as  described  above,  is  considered  the  healthy  (undamaged) 
configuration  in  our  study.  Two  problems  based  on  this  structure  are  presented  to  illustrate 
the  characteristics  of  the  procedure  of  applying  the  MRPT  to  simultaneously  determine  the 
damage  extent  in  all  property  matrices  of  undamped  (M,  K)  and  proportioned  damped  (M,  D, 
K)  structures. 


m  m 


EA  =  1.525x10''  psi        m  =  2.3  Ibs-sec^/in        Mi  =  M2  =  47.4  Ibs-sec^/in 


Figure  7.38      The  Eight-Bay  Two-Dimensional  Mass-Loaded 
Cantilevered  Truss. 


7.7.2     Proportionally  Damped  Configuration:  Damage  of  Small  Order  of  Magnitude 


In  this  first  problem,  it  is  assumed  that  the  damping  of  the  truss  is  proportional  and  equal 
to  lxlO~'  times  the  "healthy"  mass  matrix  plus  1x10"^  times  the  "healthy"  stiffness  matrix. 
In  this  example  the  damage  is  simulated  by  a  10%  stiffness  reduction  of  the  darkened  strut 
and  a  5%  reduction  of  mass  Mi.  The  damping  of  the  "damaged"  model  is  also  assumed 
proportional  and  equal  to  lxlO~'  times  the  "damaged"  mass  matrix  plus  1x10"^  times  the 
"damaged"  stiffness  matrix.  For  the  damage  analysis,  it  is  assumed  that  only  the  first  four 
modes  of  vibration  are  available.  Each  available  mode  consists  of  an  eigenvalue  and  an 
eigenvector  with  entries  at  all  FEM  DOFs.  This  present  problem  is  similar  to  a  problem 
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investigated  by  Ricles  (1991)  with  the  exception  that  in  Ricles  (1991)  the  truss  was  assumed 
undamped. 

7.7.2.1    Damage  Location 

The  first  step  in  the  analysis  is  to  use  the  damage  location  algorithm  as  discussed  in 
Chapter  5  to  determine  the  location  of  damage.  The  left  plot  of  Figure  7.39  represents  the 
unit  normalized  exact  cumulative  damage  location  vector  (CDLV)  (Eqs.  (5.2.4b)  and 
(5.4.2))  computed  from  the  exact  perturbation  matrices  AMj,  ADd  and  AIQi.  The  unit 
normalized  CDLV  computed  from  the  subspace  damage  location  algorithm  (Eqs.  (5.2.4b) 
and  (5.4.2))  is  shown  in  the  right  plot  of  Figure  7.39.  As  previously  proven,  the  location 
algorithm  is  able  to  exactly  locate  the  damage  when  presented  with  noise  free  data. 
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Figure  7.39      Problem  7.7:  Cumulative  Damage  Location  Vector:  First  Four  Modes. 


7.7.2.2    Decomposition  of  Matrix  B 

The  next  step  is  to  decompose  matrix  B  as  defined  by  Eq.  (6.4.28)  into  Bm,  Bj,  and  B^. 
As  discussed  in  Section  6.4.4,  the  proposed  decomposition  requires  that  matrix  B  be  of  full 
rank.  A  simple  singular  value  decomposition  of  matrix  B  computed  using  all  4  available 
modes  shows  that  its  rank  is  3.  Therefore,  only  three  of  the  four  available  modes  should  be 
used.  A  further  investigation  reveals  that  any  three  of  the  four  modes  result  in  a  full  rank  B. 
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Hence,  any  combination  of  three  modes  can  be  used  in  the  decomposition.  The  three  plots  in 
the  upper  half  of  Figure  7.40  show  the  unit  normahzed  exact  cumulative  vector  (Eq.  (6.4.35), 
Section  6.4.4)  computed  using  the  exact  perturbation  matrices,  AMj,  ADj  and  AK^,  and  the 
first  three  modes.  The  unit  normalized  cumulative  vector  associated  with  the  computed  Bm, 
Bd,  and  B^  using  the  proposed  decomposition  procedure,  for  the  first  three  modes,  are 
displayed  in  the  lower  half  of  Figure  7.40.  A  comparison  between  upper  and  lower  plots 
shows  that  the  decomposition  procedure  resulted  in  the  exact  Bm,  Bj,  and  Bk  for  the 
associated  modes.  This  is  also  true  if  any  other  combination  of  three  modes  was  used.  Note 
that  the  cumulative  vectors  in  Figure  7.40  give  information  on  the  effect  of  the  damage  on  the 
mass,  damping  and  stiffness  properties. 
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Figure  7.40      Problem  7.7:  Cumulative  Vectors  Associated  With  the 
Exact  and  Computed  Bm,  Bj,  B^:  First  Three  Modes. 
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7.7.2.3    Damage  Extent 

Now  that  Bm,  Bd,  and  B^  have  been  computed,  the  final  step  is  to  compute  the 
perturbation  to  the  property  matrices  due  to  the  damage  (AMj,  ADj,  AK^)  using  the 
minimum  rank  perturbation  theory  (MRPT)  as  defined    in  Section  6.4.4.    Using  the 
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Figure  7.41      Problem  7.7:  Exact  and  Computed  AMd,  ADd,  AIQ 


connectivity  information  of  the  original  FEM  along  with  the  information  of  the  damage 
effect  on  the  three  property  matrices  (shown  in  Figure  (7.39)),  it  can  be  deduced  that  the  rank 
of  the  "true"  AMd  is  two  and  the  rank  of  the  true  AKd  is  one.  Because  of  the  fact  that  the 
damage  in  the  damping  is  proportional  to  the  damage  in  the  mass  and  stiffness,  it  can  also  be 
deduced  that  the  rank  of  the  "true"  ADd  is  three.  The  rank  information  is  ascertained  by 
identifying  which  structural  elements  connect  the  damaged  DOFs,  and  then  adding  up  the 
rank  of  each  damaged  structural  element.  The  rank  estimation  could  also  be  obtained  by 
performing  a  singular  value  decomposition  on  B^,  Bd  and  B^.  From  the  rank  information 


153 


and  propositions  6. 1  and  6.2,  it  is  clear  that  only  two  modes  are  needed  to  computed  AMd 
using  the  MRPT.  Likewise,  only  three  modes  are  needed  to  computed  ADj  and  only  one 
mode  is  needed  to  compute  AK^.  Since  we  are  dealing  with  noise  free  eigendata,  any  two, 
three  modes,  or  any  one  mode  can  be  used  in  the  computation  AMd,  ^Dd  or  AKd, 
respectively.  Mesh  plots  of  the  exact  and  computed  AMd,  ^Dd  and  AIQ  matrices  are  shown 
in  Figures  7.41.  Again,  a  comparison  between  exact  and  computed  AMd,  ^Dd  and  AKd 
shows  that  the  MRPT  is  able  to  reproduce  the  exact  damage  effect. 

In  this  problem,  as  shown  in  this  investigation,  the  proposed  damage  detection 
technique  is  able  to  assess  the  damage  exactly  when  three  noise  free 
eigenvalues/eigenvectors  are  used.  Unfortunately,  in  practice,  the  measured  eigendata  are 
always  corrupted  by  noise.  To  simulate  some  kind  of  a  practical  situation,  random  noise  was 
added  to  the  eigendata.  When  using  noisy  eigendata  (  first  four  modes),  the  present 
algorithm  was  unable  to  locate  the  damage.  The  reason  is  that  the  damage  as  simulated  is  of 
small  order  of  magnitude  which  results  in  only  a  small  change  in  the  eigendata  of  the  healthy 
model.  In  this  case,  the  noise  totally  masks  the  damage  as  reflected  in  the  "measured"  noisy 
eigendata. 

7.7.3     Undamped  Configuration:  Damage  of  Large  Order  of  Magnitude 

In  this  problem,  the  truss  is  assumed  undamped.  The  damage  is  simulated  by  a  60% 
stiffness  reduction  of  the  darkened  strut  and  a  40%  reduction  of  mass  Mi .  For  the  analysis,  it 
is  assumed  that  only  the  first  five  modes  of  vibration  are  available.  Each  available  mode 
consist  of  an  eigenvalue  and  an  eigenvector  with  entries  at  all  FEM  DOFs.  This  problem  is 
investigated  for  four  different  cases.  Each  case  represents  a  different  amount  of  random 
noise  added  to  the  eigenvectors  of  the  available  modes:  0  %  ,  2.5  %,  5%  and  10%. 

7.7.3.1    Noise  Free  Eigendata 

The  procedures  used  in  the  investigation  of  this  case  are  similar  to  the  one  used  in  the 
previous  problem.     The  computed/exact  cumulative  damage  location  vector  and  the 


154 


cumulative  vector  associated  with  Bm  and  B^  are  displayed  in  the  first  row  of  plots  of  Figure 
7.42.  Likewise,  the  first  row  of  plots  of  Figure  7.43  show  the  computed/exact  AMj  and  AKj 
when  using  noise  free  data.  Again,  as  previously  shown,  exact  data  provides  exact  results. 

7.7.3.2    Noisy  Eigendata 

The  effect  of  introducing  noise  into  the  "measured"  eigenvectors  is  shown  in  the  lower 
three  rows  of  plots  in  Figures  7.42  and  7.43.  In  Figure  7.42,  it  is  obvious  that  as  noise  is 
increased,  the  cumulative  damage  vector  B  (first  column)  becomes  corrupted,  to  the  point 
that  when  10%  noise  is  added  it  is  difficult  to  ascertain  information  concerning  the  state  of 
damage.  In  performing  the  decomposition  of  B  into  B^,  and  B^  one  should  only  use  those 
modes  of  vibrations  whose  damage  vector  dj  reflects  the  same  "nature  of  damage"  as  the 
cumulative  damage  vector  B.  For  this  particular  example,  only  modes  3, 4  and  5  meet  this 
criterion  for  the  cases  of  2.5%  and  5%  random  noise.  For  the  1 0%  noise  case,  the  cumulative 
damage  vector  and  the  individual  mode  damage  vectors  did  not  clearly  indicate  damage.  To 
maintain  a  level  of  consistency,  the  same  mode  set  (3,  4,  5)  was  used  for  the  case  of  10% 
noise.  The  decomposition  of  B  into  Bm  and  B^  for  the  case  of  noisy  measurements  is  shown 
in  the  remaining  rows  and  columns  of  Figure  7.42.  As  shown  in  Figure  7.42,  the  added  noise 
has  two  effects  in  the  decomposition  process.  The  first,  which  is  rather  obvious,  is  that  noise 
causes  the  appearance  of  small  numerical  component  at  all  DOFs,  although  in  the  low  level 
noise  cases  (2.5%  and  5%)  it  is  still  clear  which  DOFs  have  been  actually  damaged.  The 
second  more  subtle  effect  is  that  the  actual  decomposition  process  is  no  longer  exact.  This 
manifests  itself  in  that  the  computed  B^  has  substantial  indicators  of  damage  at  DOFs  in 
which  there  is  actually  only  stiffness  damage.  The  same  is  true  for  Bk,  in  that  damage  is  also 
indicated  at  DOFs  in  which  there  is  only  mass  damage.  These  indicators  are  shown  darkened 
in  Figure  7.42. 

Figure  7.43  shows  the  effect  of  noise  on  the  computed  perturbation  matrices.  Mode 
selection  was  based  on  (i)  rank  calculation  of  Bm  and  B^,  and  (ii)  comparison  of  individual 
columns  of  Bm  and  Bk  to  their  associated  cumulative  vector.  For  the  cases  of  2.5%  and  5% 
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noise,  the  extent  calculation  stills  provides  a  reasonable  indication  of  damage,  although  the 
magnitude  of  damage  is  in  error.  The  percentage  errors  of  the  estimated  mass  and  stiffness 
damage  with  respect  to  the  exact  damage  are  listed  in  Table  7.11.  Figure  7.43  shows  the 
extent  calculation  with  no  enhancements. 

For  the  case  in  which  10%  noise  was  added,  the  damage  location  and  extent  algorithm 
failed.  In  actual  practice,  the  analysis  would  have  halted  after  the  calculation  of  the 
individual  and  cumulative  damage  vectors.  At  that  point,  inspection  of  the  damage  vectors 
would  give  no  clear  indication  as  to  the  state  of  damage.  The  decomposition  of  B  into  Bm  and 
Bic  as  well  as  the  calculation  of  the  perturbation  matrices  was  just  carried  out  for 
completeness  of  the  results. 


Table  7. 1 1        Problem  7.7:  Percentage  Error  of  Damage  Estimate  with 
respect  to  Exact  Damage 


Percentage  Error  With  Respect  to 
Exact  Damage 

Noise  Added 

Mass 

Stiffness 

0.0% 

0.00 

0.00 

2.5% 

6.60 

9.10 

5.0% 

99.10 

24.60 

7.8    Summary 


A  technique  that  approaches  the  structural  damage  assessment  problem  in  a  decoupled 
fashion  was  demonstrated  and  evaluated  using  numerical  and  actual  experimental  test  data. 
In  all  presented  examples,  the  structural  damage  was  first  located  by  using  the  subspace 
rotation  damage  location  algorithm  formulated  in  Chapter  5.  With  location  determined,  the 
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minimum  rank  perturbation  theory  (MRPT)  developed  in  Chapter  6  is  then  applied  to  assess 
the  extent  of  the  damage.  As  illustrated,  the  decomposition  of  the  damage  assessment 
problem  has  two  distinct  advantages.  First,  each  subproblem  is  shown  to  be  computationally 
attractive.  Second,  the  decoupling  allows  engineering  judgement  to  enter  into  the  extent 
algorithm.  By  making  use  of  the  results  of  the  damage  location  algorithm,  the  modes  which 
do  not  exhibit  the  damage  may  be  eliminated  before  the  application  of  the  MRPT  extent 
algorithm.  In  general,  the  subspace  rotation  damage  location  algorithm  and  the  MRPT 
performed  well  in  assessing  damage  in  all  studied  examples.  The  eigenvector  filtering 
algorithm  developed  in  Chapter  5  was  shown  to  be  useful  in  improving  the  accuracy  of  the 
damage  extent  assessment.  The  MRPT  was  also  demonstrated  to  be  applicable  in  model 
refinement  problems. 


CHAPTER  8 
CONCLUSION  AND  SUGGESTIONS  FOR  FUTURE  WORK 


This  study  investigated  the  development  of  four  algorithms  relevant  to  the  areas  of  finite 
element  model  refinement  and  structural  damage  assessment.  All  four  algorithms  make  use 
of  an  original  finite  element  model  (FEM)  and  a  subset  of  experimentally  measured 
eigenvalues  and  eigenvectors.  The  first  algorithm,  termed  the  inverse/hybrid  approach,  is  a 
model  refinement  algorithm  with  basis  in  the  standard  framework  of  the  inverse  problem.  In 
the  formulation  of  the  inverse/hybrid  approach,  it  is  proposed  to  approximate  the 
unmeasured  experimental  modes  with  the  corresponding  analytical  modal  information. 
With  a  complete  hybrid  set  of  experimental  and  analytical  modal  parameters,  the  symmetry 
of  the  inaccurate  finite  element  property  matrices  is  enforced  by  mass  orthogonalizing  the 
eigenvectors.  An  orthogonalization  strategy  that  assigns  more  credibility  to  the  measured 
eigenvectors  was  discussed.  Through  an  example,  it  was  shown  that  the  performance  of  the 
inverse/hybrid  approach  is  similar  to  the  algorithm  formulated  by  Baruch  and  Bar  Itzhack 
(1978). 

Next,  a  subspace  rotation  technique  to  improve  the  performance  of  an  existing  model 
refinement  algorithm,  termed  SEAMRA  (Zimmerman  and  Widengren,  1989),  was 
presented.  The  mathematical  formulation  of  SEAMRA  is  based  on  a  control  concept  known 
as  the  eigenstructure  assignment.  The  proposed  improvement  results  in  both  an 
enhancement  of  the  eigenvector  assignability  and  a  substantial  decrease  of  the  computational 
burden.  The  enhanced  SEAMRA  was  shown  to  be  suited  for  both  model  refinement  and 
structural  damage  assessment. 

An  efficient  structural  damage  location  algorithm  that  bypasses  the  general  framework 
of  the  model  refinement  problem  was  then  formulated.   The  algorithm,  also  termed  the 
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subspace  rotation,  is  similar  to  the  Modal  Force  Error  Criteria  formulated  by  several 
researchers.  The  effects  of  measurement  error  on  the  performance  of  the  subspace  rotation 
were  discussed.  A  new  viewpoint  that  reduces  the  effects  of  measurement  error  for  certain 
classes  of  structure  was  developed.  Additionally,  the  use  of  multi-mode  measurements  as  a 
technique  to  overcome  the  effects  of  these  errors  is  also  discussed.  Furthermore,  an 
eigenvector  filtering  technique  was  formulated. 

A  minimum  rank  perturbation  theory  (MRPT)  was  then  developed.  A  number  of 
computationally  attractive  damage  extent  algorithms  based  on  the  MRPT  were  formulated. 
The  formulations  of  the  presented  MRPT  based  algorithms  are  consistent  with  the  effects  of 
most  structural  damage  on  FEMs.  The  formulated  algorithms  address  several  damage 
scenarios  for  undamped,  proportionally  damped  and  nonproportionally  damped  FEMs.  A 
damage  scenario  is  defined  by  the  effect  of  the  structural  damage  on  the  FEM.  Two 
techniques  to  improve  the  performances  of  the  MRPT  based  extent  algorithm  were 
discussed.  The  first  technique  suggests  the  exclusion  in  the  damage  extent  calculations  of  the 
measured  "damaged"  modes  that  do  not  reflect  the  state  of  the  damage  as  determined  by  the 
subspace  rotation  damage  location  algorithm.  Alternatively,  the  eigenvector  filtering 
algorithm  discussed  earlier  was  shown  to  be  also  useful  in  improving  the  damage  estimate. 

Finally,  the  subspace  rotation  damage  location  algorithm  and  the  MRPT  based  damage 
extent  algorithms  were  evaluated  using  both  computer  simulated  and  actual  experimental 
data.  All  issues  raised  in  the  formulation  of  these  algorithms  were  demonstrated.  In  every 
problem,  the  damage  assessment  process  was  approached  in  a  decoupled  fashion.  The 
location  of  the  damage  was  first  determined.  The  extent  of  the  damage  was  then  assessed  by 
making  use  of  the  results  of  the  damage  location.  The  decoupled  approach  improves  both  the 
accuracy  and  the  efficiency  of  the  extent  computations.  In  general,  all  evaluated  algorithms 
performed  very  well  in  assessing  the  damage  in  the  examples  studied.  The  algorithms  show 
great  promises  in  handling  "real  life"  structures.  This  was  demonstrate  in  the  investigation  of 
the  experimental  data  and,  in  particular,  the  NASA  8-bay  truss  examples. 
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To  fully  demonstrate  the  practicality  of  the  MRPT  based  damage  extent  algorithms  and 
the  subspace  rotation  damage  location  algorithm  in  handling  "real  life"  structures,  the  issue 
of  incomplete  eigenvector  measurement  must  be  intensively  investigated.  This  problem  was 
pointed  out  in  Chapter  2.  Two  alternative  approaches,  model  reduction  and  eigenvector 
expansion,  were  discussed  as  possible  solutions  to  the  incomplete  measurement  problem. 
These  approaches  have  only  been  investigated  in  the  context  of  the  model  refinement 
problem.  An  investigation  of  their  performance  vis-a-vis  the  damage  assessment  problem  is 
of  order.  Additionally,  formulations  of  damage  assessment  algorithms  that  do  not  require  the 
use  of  a  finite  element  model  should  be  considered.  This  could  be  useful  when  dealing  with 
"old"  structures  that  have  no  analytical  model. 
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